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ABSTRACT 

This paper presents the measurement of the neutron star (NS) radius using the thermal spectra 
from quiescent low- mass X-ray binaries (qLMXBs) inside globular clusters (GCs). Recent observations 
of NSs have presented evidence that cold ultra dense matter - present in the core of NSs - is best 
described by "normal matter" equations of state (EoSs). Such EoSs predict that the radii of NSs, 
i?NSi are quasi-constant (within measurement errors, of ~ 10%) for astrophysically relevant masses 
(Mns > 0.5 Mq). The present work adopts this theoretical prediction as an assumption, and uses it 
to constrain a single i?NS value from five qLMXB targets with available high signal-to-noise X-ray 
spectroscopic data. Employing a Markov-Chain Monte-Carlo approach, we produce the marginalized 
posterior distribution for j?NS> constrained to be the same value for all five NSs in the sample. An 
effort was made to include all quantifiable sources of uncertainty into the uncertainty of the quoted 
radius measurement. These include the uncertainties in the distances to the GCs, the uncertainties 
due to the Galactic absorption in the direction of the GCs, and the possibility of a hard power-law 
spectral component for count excesses at high photon energy, which are observed in some qLMXBs 
in the Galactic plane. Using conservative assumptions, we found that the radius, common to the 
five qLMXBs and constant for a wide range of masses, lies in the low range of possible NS radii, 
-Rns = 9.lj^'Jkm (90%-confidence). Such a value is consistent with low-i?NS equations of state. We 
compare this result with previous radius measurements of NSs from various analyses of different types 
of systems. In addition, we compare the spectral analyses of individual qLMXBs to previous works. 
Subject headings: stars: neutron — X-rays: binaries — globular clusters: individual (w Cen, M13, 
M28, NGC 6397, NGC 6304) 



1. INTRODUCTION 

The relation between pressure and energy den- 
sity in matter at and above the nuclear saturation 
density p r . = 2.8 x 10 14 gem" 3 is largely unknown 



(jLattimer fc Prakashl l200ll l2007h . This is mostly due 
to uncertainties of many-body interactions as well as 
the unknown nature of strong interactions and symme- 
try energy. Inside neutron stars (NSs), the equation of 
state of dense matter (P(e), written dEoS, hereafter) 
can be mapped into a mass-radius relation A/ns(-Rns) 
by solving the Tolman-Oppenheimer-Volkoff eq uation 
(|Oppenheimer fc Volkoffl[l939l: iMisner et al.llT973l ). His- 
torically, well before any observational constraints could 
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be placed on the dEoS, nuclear theory attempted to de- 
termine the P (e) relation that would govern the behav- 
ior of cold ultra-dense matter. Since the cores of NSs are 
composed of such matter, its behavior is of astrophysical 
interest; likewise, the behavior of NSs due to the compo- 
sition of its core is of nuclear physics interest. 

Three main families of dEoSs have been discussed 
in the last 10-20 years. The first one regroups "nor- 
mal" dense matter EoSs. At densities at p c , nuclei dis- 
solve and merge, leaving undifferentiated nuclear mat- 
ter in /3-equilibrium. In this type of matter, the pres- 
sure is neutron-dominated via the strong force, with a 
small proton fraction. In other words, NSs are pres- 
sure supported against gravity by neutron degeneracy. 
The "normal" dEoSs are calculated with a relativis- 
tic treatment of nucleon-nucleon interactions, leading 
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to a relation between pressure a nd density, with the 
press ure vanishing at zero densities (|Lattimer fc Prakashl 
I2001f) . For NSs, such dEoSs correspond to M NS (i?Ns) 
lines composed of two parts. One corresponds to con- 
stant low A/ns at large i?NS values. Then, as the den- 
sity increases, the Mns(-Rns) relation for "normal" dE- 
oSs evolves to quasi-constant 1 i?NS as Mns increases, 
up to a maximum mass, above which the NS collapses 
to a black-hole. Examples of the proposed form of 
these dEoSs i nclude AP3-AP4 (lAkmal fe Pandharipande] 
1997ft . ENG (lEngvik et alj|1996ft. MP A l (iMuther et alj 
1981. MSO and MS2 (IMuller fe Serotl 11996ft. and LS 
(jLattimer fe Swestvlll991ft . 

A second family of dEoSs is characterized by matter 
in which a significant amount of softening (i.e., less pres- 
sure) is included at high densities, due usually to a phase 
transition at a critical density which introduces an ad- 
ditional hadronic or pure-quark component in what is 
referred to as the NS's "inner core" . Additional compo- 
nents, such as a population of hyperons at large densities 
(GM3, IGlendenning fe M oszkowski 1991) , or kaon con- 
densa tes (GS1, GS2, Glcndcnning & Schaffncr- Bielichl 
1999), have been considered. For that reason, these dE- 
oSs are referred to as "hybrid" dense matter. Because 
of this phase-transition, the maximum Mns is rather 
low (Mns < 1-7 M Q ). This also results in Mns(-Rns) 
curves with a smooth decrease in Mns from the maxi- 
mum to the minimum Mns , as /?ns increases. So me of 
the "hybrid'' dEoS s are MSI (IMuller fc SerotllT996ft . FSU 
(IShen et al.l l 20l0allbT). GM3 (IGlendenning fc Moszkowski 
11991ft . GS1 (IGlendenning fc Schaffner-Bielichlll999h . and 
PAL6 (jPrakash et al.lll988ft . 

The third family of dEoSs relies on the assumption that 
strange quarks compose matter in its ground state. One 
characteristic of such matter is that the pressure van- 
ishes at a non-zero density, compared to the other types 
of matter described above - they have solid surfaces. In 
Mns-Rns space, these quark star dEoSs follow lines of 
increasing Mns with increasing radius, up to a maxi- 
mum radius. Above this value, i?NS starts decreasing 
as Mns increases until Mns reaches its own maximum, 
where the object collapses to a black- hole. The maxi- 
mum i?NS varies between ~ 9 km and ~ 11km, depend- 
ing on the model parameters used, namely, the strange 
quark m ass m s and the quan tum chromodynamic cou- 
pling a c (jPrakash et alj 11995ft . Note that "hybrid" and 
"normal" matter stars do not have this constraint, and 
their radii can theoreti cally be as large as ~ 100 k m, at 
masses M NS < 0.5 M (jLattimer fc Prakash|[200TI) 

Since matter at such densities cannot be produced in 
Earth laboratories, constraints on the dEoS theoretical 
models can only be placed by the study of NSs, the 
only objects in the Universe containing matter at such 
densities. The measurements of Mns and i?NS have 
the potential to provide great insight to the theory of 
cold ultra dense matter. Various methods exist to m ea- 
sure Mns and i?NS (e.g.. ILattimer fc Pra kash 20Q3, for 
a general review). These include the stu dy of quasi- 
periodic oscillations in active X-ray binaries (jMiller et al.1 

1 Here, and elsewhere, we use the term "quasi-constant" to mean 
constant within measurement precision, ~ 10%. This should be 
differentiated from a value which is constant when measured with 
infinite precision, or a value which is constant according to theory. 



119981: iMendez fc Bellonil 12007ft. Keplerian parameters 
in NS binaries (|Nice et al.1 12004 iDemorest et~aTI 120101 

for A/ns measurements), thermonuclear X-ray bursts 
(|Ozelll2006t ISuleimanov et al.ll2011bl for M NS -.%s mea- 
surements), pulse-timing analysis of m illisecond pulsars 
(jBogdanov et al.ll2008l: [Bogdanovll2012ft . and the thermal 
spectra of quiescent low-mass X-ray binaries (qLMXBs), 
which is the method of this investigation. Each of these 
different methods have their own unique systematic un- 
certainties, and it is therefore of value to pursue each, to 
permit intercomparison of their conclusions. 

By itself, a A/ns measurement can only place new con- 
straints on the dEoS when the measured value is above 
that of all previous Mns measurements. In A/ns~-Rns 
space, each dEoS predicts a maximum Mns, above which 
the NS collapses to a black hole. In particular, hybrid 
dEoSs are characteriz ed by a relatively low ma ximum 
M NS (Mns < 1.8M ? . lLa"ttimer fc Prakashll2001ft . while 
normal matter dEoSs produce ma ximum Mns of up to 
2.5 M© (jLattimer fc Pra kash 200j|). The maximum M N s 
for strange quark matter (SQM) dEoSs is typ ically in the 
vicinity of 2M (jLattimer fc Prakashl [20011 ) . The max- 
imum Mns property of EoSs can be used to excludes 
dEoSs. Historically, Msrs measurements were in the 1.3- 
1.5 Mq range. While the first precise Mns measure- 
ments confirmed theoreti cal predictions about NSs (e.g., 
iTavlor fc Weisberg] 11989ft . subsequent measurements at 
and below previous values did not place any new con- 
straints on the dEoS. Recently, the mass of the radio pul- 
sar PSR 1614—2230 was precisely measured w ith a value 
Af N s = 1-97 ± 0.04 M Q (jDemorest et al.|[2010ft . The im- 
plications of this measurement f or nuclear physi cs have 
been discussed with some depth iLattimerl (|2011ft . Such 
a high Mns excludes previously published hybrid mod- 
els of dEoSs (using specific values of assumed parameters 
from within their allowed regions), although it does not 
rule out any specific form of exotica. SQM dEoSs also 
seem to be disfavored, since their predicted maximum 
Mns approaches the 2 M Q limit for only some of the 
models within the parameter spaces permitted by nuclear 
physics constraints. Nonetheless, fine tuning of models 
may allow these disfavored dEoSs to be marginally con- 
sistent w ith the M NS measured in PSR 1614-2230 (for 
example, iBednarek et al.ll201lUWeissenborn et al.ll2012l 
for hybrid models. and lLai fc Xull201lL for SQM models). 

Overall, this high-A/NS measurement seems to favor 
"normal matter" hadronic dEoSs. This would mean that 
the radius of astrophysical NSs should be observed to 
be within a narrow ( < 10%) range of values for Mns > 
0.5 Mq, since "normal matter" dEoSs follow lines of 
quasi-constant radius in AfN S--RNS-space at such masses 
(Lattimer & Prakash 2001). It is important to notice 
that the spread in i?NS increases for stiff EoSs, espe- 
cially close to the maximum A/ns of the compact object 
(e.g., up to a 2-km difference in i?NS for the EoS PALI, 
iPrakash et aTl[l9"88ft . 

The empirical dEoS obtained from Mns - -Rns confi- 
dence regions from type-I X-ray bursts and from the 
thermal sp ectra of qLMXBs combin ed also favors this 
conclusion (jSteiner et al.l 120101 12012). Using a Bayesian 
approach, the most probable dEoS was calculated, re- 
sulting in a dEoS approaching the behavior of theoret- 
ical hadronic dEoSs, with predicted radii in the range 
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i?NS ~ 10 — 13 km. Such radii suggest that soft hadronic 
dEoSs are describing the dense matter inside NSs. How- 
ever, different analyses of other NSs found radii consis- 
tent with stiff dEoSs. These include the qLMXB X7 in 
47Tuc (IHeinke et all [20061) . or the type I X -ray burster 
4U 1724-307 ([Suleimanov et all 1201 lal ). Nonetheless, 
these results are not inconsistent with the observation 
that the i?NS is almost constant for a large range of -Mns, 
since they are consis tent with stiff "norma l matter" dE- 
oSs, such as MSO/2 (IMtiller fe Serotl[l99"l 

Given the evidence supporting the "normal matter" 
hadronic dEoSs, it therefore becomes a natural assump- 
tion - to be tested against data - that observed NSs 
have radii which occupy only a small range of i?NS val- 
ues ( $ 10%). Using the thermal spectra of five qLMXBs, 
fitted with a H-atmosphere model, a single _/?ns value is 
assumed and measured, as well as its uncertainty. Fur- 
thermore, under this assumption, the best-fit Mns and 
surface effective temperature kT c g for these qLMXBs and 
their uncertainties are extracted. The various sources 
of uncertainty involved in this spectral analysis are ad- 
dressed, including, the distances to the qLMXBs, the 
amount of galactic absorption in their direction, and the 
possibility of an excess of high-energy photons as ob- 
served for other qLMXBs (and modeled with a power-law 
component, PL hereafter). The goal is to place the best 
possible constraints on _Rns accounting for all know un- 
certainties, and eliminating all unquantifiable systematic 
uncertainties. 

In this article, we provide the necessary theoretical 
background and observational scenario to understand 
Mns - -Rns measurements of NSs from qLMXBs (§ [5]). 
The organization of the rest of the paper is as follows: 
Section [3] explains the analysis of the X-ray data. Sec- 
tion @] contains the results of the spectral analysis. A 
discussion of the results is in Section [5] and a summary 
is provided in Section [5] 



2. QUIESCENT LOW-MASS X-RAY BINARIES 

The low-luminosity of qLMXBs was initially observed 
following the outbursts of the X-ray t ransients Cen X-4 
and Aql X-l ([van Paradiis e t al. 1983). This faint emis- 
sion (Lx ~ 10 32 — 10 33 erg s" 1 , 4-5 orders of magnitude 
fainter than during outburst) was originally interpreted 
as a thermal blackbody emission. Low-level mass accre- 
tion onto the compac t object was thought to explain the 
observed luminosity ([Verbunt et al.lll994l ). 

Later, an alternative to the low-level accretion hypoth- 
esis was proposed. This alternate theory, which be- 
came the dominant explanation for the emission from 
qLMXBs, suggests that the observed luminosity is pro- 
vided, not by low-M, but by t he heat deposited in the 
deep crust during outbursts ([Brown et al.l [1998). In 
the theory of deep crustal heating (DCH), the mat- 
ter accreted during an outburst releases ~ 1.9 McV 
of energy via pressure-sensitive reactions: electron 
captures, neutron emissions or pycnonuclear reactions 
<tSatolH979t IHaensel fe Zduifil [19901 IGupta et al J 120071: 
Hae nsel fc Zdunik 2008). Therefore, the time-averaged 
quiescent luminosity is proportional to the time-averaged 



mass accretion rate: 



(L) = 9x10 



32 



(M) 



Q 



10 



M Q yr- 1 1.5MeV/amu 



erg s 



(1) 

where Q is the avera ge heat deposited i n the NS cru st 
per accreted nucleon ([Brown et al.l fl998: Bro wnll2000[ ). 

Following this hypothesis about the energy source of 
the quiescent luminosity, the theory of DCH also ex- 
plains the observed spectra of qLMXBs. As a result of 
the energy deposited in the deep crust, the core heats up 
during the outbursts. The energy is then re-radiated 
on core-cooling time scales a way from the crus t, and 
through the NS atmosphere ([Brown et alJ Il998f >. The 
NS atmosphere is assumed to be composed of pure hy- 
drogen. Indeed, at the accretion rates expected during 
quiescence, he avier elements settle on time scales of or- 
der ~ seconds (|Romani|11987t iBildsten et al.l ll992). The 
possibility of helium (He) or carbon atmospheres around 
NSs in LMXBs has a lso been studied (|Ho fc Heinkell2009t 
IServillat et aII[20Tl . 

Several m odels of H-atmosphere a round N Ss have been 
devel op ed dRaiagopal fc Romaml 119961 iZavlin et all 
[1991 iMcClintock et all 120041: IHeinke et all 120061: 
IHaakonsen et alll2012l) . Thev are now routinely used to 
explain the emergent spectra of qLMXBs, with emission 
area radii compatible with the entire surface area of 
NSs, compared to derived emissi on area radii of < 1 km 
in the blackbody interpretation (|Rutledge et al.l fl999[) . 
The DCH theory and H-atmosphere models were first 
applied to explain the quiescent spectra and measure the 
radius of historically transient LMXB (e.g., Cen X-4, 
Campana efaTI 120001 IRutledge etHI I2001al Aql X-l, 



Rutledge et alj|2001b ). However, the 10-50% systematic 



uncertainty on the distance to field LMXBs directly 
contributes to a 10-50% uncertainty on the radius mea- 
surements. Due to these large systematic uncertainties, 
these objects provide limited use to place constraints 
on the dEoS until more precise measurements of their 
distances can be obtained. 

Placing tight constraints on the dEoS requires ~ 5% 
uncertainty on the -Rns measurements. This constraint 
is approximately the half-width of the range of radii in 
the Mns (-Rns) relationships corresponding to "normal 
matter" EoSs. Globular clusters (GCs) have proper- 
ties which make them ideal targets for qLMXB obser- 
vations: relatively accurately measured distances; bet- 
ter characterized Galactic absorption; over-abundances 
of LMXBs; and LMXBs with magnetic field weak enough 
that the thermal spectrum is not affected ([Heinke et all 
2006). A handful of qLMXBs have been discovered in 
GCs so far; only a few have X-ray spectra with the high 
signal-to- noise ratio (S/N) necessary to measure i?NS 
with < 10 — 15% uncertainty, including the uncertainty 
to their distances. 

3. DATA REDUCTION AND ANALYSIS 

3.1. Targets 

The targets used in this work are chosen among the 
qLMXBs located in GCs that produced the best i?NS 
measurements, i.e., with i?oo uncertainties of ^5 15% in 
the previous works. 

The GCs w Ce n (IRutledge et all [20021: iGendre et al 
I2003al) and M13 (|Gendre et al.ll2003bl : iCatuneanu et al 
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TABLE 1 

X-ray Exposures of the Targeted Clusters. 



Target 


Obs. ID 


Starting 


Usable time 


S/N 


Telescope 


Filter 


Refs. 






Time fTT^ 

i line ^ r l j 


(kscc) 




and detector 


or AJode 




M28 


2683 


2002 July 04 18:02:19 


14.0 


23.85 


Chandra ACIS-S3 (BI) 


VFAINT 


2 


M28 


2684 


2002 Aug. 04 23:46:25 


13.9 


23.54 


Chandra ACIS-S3 (BI) 


VFAINT 


2 


M28 


2685 


2002 Sep. 09 16:55:03 


14.3 


23.90 


Chandra ACIS-S3 (BI) 


VFAINT 


2 


M28 


9132 


2008 Aug. 07 20:45:43 


144.4 


78.75 


Chandra ACIS-S3 (BI) 


VFAINT 


1, 3 


M28 


9133 


2008 Aug. 10 23:50:24 


55.2 


48.46 


Chandra ACIS-S3 (BI) 


VFAINT 


1, 3 


NGC 6397 


79 


2000 July 31 15:31:33 


48.34 


25.03 


Chandra ACIS-I3 (FI) 


FAINT 


4, 5 


NGC 6397 


2668 


2002 May 13 19:17:40 


28.10 


25.47 


Chandra ACIS-S3 (BI) 


FAINT 


5 


NGC 6397 


2669 


2002 May 15 18:53:27 


26.66 


24.97 


Chandra ACIS-S3 (BI) 


FAINT 


5 


NGC 6397 


7460 


2007 July 16 06:21:36 


149.61 


52.31 


Chandra ACIS-S3 (BI) 


VFAINT 


5 


NGC 6397 


7461 


2007 June 22 21:44:15 


87.87 


41.40 


Chandra ACIS-S3 (BI) 


VFAINT 


5 


M13 


0085280301 


2002 Jan. 28 01:52:41 


18.8 


14.25 


XMM pn, MOSl, M0S2 


Medium 


6,7 


M13 


0085280801 


2002 Jan. 30 02:21:33 


17.2 


12.10 


XMM pn, MOSl, MOS2 


Medium 


6,7 


M13 


5436 


2006 Mar. 11 06:19:34 


27.1 


16.07 


Chandra ACIS-S3 (BI) 


FAINT 


1,8 


M13 


7290 


2006 Mar. 09 23:01:13 


28.2 


16.01 


Chandra ACIS-S3 (BI) 


FAINT 


1,8 


uj Cen 


653 


2000 Jan. 24 02:13:28 


25.3 


13.33 


Chandra ACIS-I3 (FI) 


VFAINT 


9 


a) Cen 


1519 


2000 Jan. 25 04:32:36 


44.1 


16.45 


Chandra ACIS-I3 (FI) 


VFAINT 


9 


u) Cen 


0112220101 


2001 Aug. 12 23:34:44 


33.9 


24.35 


XMM pn, MOSl, MOS2 


Medium 


7,10 


NGC 6304 


11074 


2000 July 31 15:31:33 


98.7 


27.94 


Chandra AC1S-13 (FI) 


VFAINT 


1 



Note. — TT refe rs to Terrestri a l Tim e. FI and BI refers to the fro nt-ill uminated and back-illu minat ed ACIS chips . References: 
(1) This work: (21 IBecker et all l|2003l); (31 IServiliaTc t al. (2012 j); (4) IGrindlav et alj J200B: (51 iGuillot et all fl2011al); (6) 
IGendre et"aH l l2003bl) : (7) IWebb &: Barrel (1200711 : (81 ICatuneanu et ail 1120131) , (91 IRutledge et all j2002T) : (10) IGendre et al.l 
; All observations have been re-processed and rc-analyzed in this work. The references provided here arc given to indicate 
the previously published analyses of the data. 




Time (see) 

Fig. 1. — This figure shows the XMM -pn full-detector light curves 
of w Cen, ObsID 0112220101, with bins of 100 sec. The black (top) 
line corresponds to the pn camera light curve. The t = sec time is 
the beginning of the pn exposure, on 2001 Aug. 12, 23:34:44. The 
MOSl light curve, in red (bottom), is shown for completeness and 
because periods of flaring are more readily visible. The time inter- 
vals with large background flaring are excluded from the analyzed 
data set. 

2013) each have one qLMXB that was used in previ- 
ous work to place m oderate constraints on the dEoS 
(jWebb fc BarretJl2007t ). The projected radius measure- 
ments Roo reported in the original works are within 2-3% 
uncertainty. However, there is evidence that these uncer- 
tainties are highly under-estimated (§2). 
The qLMXB in the core of NGC 6304, discovered re- 



cently with the XMM-Newton observatory (XMM, here- 
after) and confir med with a short Chan dra X-ray Obser- 
vatory exposure (jGuillot et al.l l2009al fbh , was then ob- 
served for 100 ks with ACIS-I onboard Chandra (Ad- 
vanced Charge-coupled-device Imaging Spectrometer). 
In this work, only the long Chandra exposure is used 
since, in the XMM observation, the core source is con- 
tamina ted by nearby source s, mostly one spectrally hard 
source (|Guilk>t et al.ll2009bf ). 

The qLMXB i n NGC 6397 (named U24 in the discov- 
ery observation, iGrindlav et al.l 120011 ) has a i?NS value 
measured with ~ 8% uncertainty, obtained from a total 
of 350 k s of Chandra X-ray Observatory archived obser- 
vations (jGuillot et al.ll201 laf) . The spectra for this target 
were re-analyzed in this work, for a more uniform anal- 
ysis. 

Finally, the Roc measurement of the NS qLMXB in 
the core of M28 reported in the discovery observation 
does not place useful constr aints on the dEoS: R^ = 
14.5±g;g (jBecker et alj 120031 ). However, an additional 
200 ks of archived observations with Chandra have been 
analyzed in a recent work, finding R^s = 9±3 km and 
Mjsss = l-C^g-M© for a H-atmosphere, and i?NS = 
14±| and M NS = 2.0 t°;^ M for a pure He-atmosphere 
(|Servillat et al.l [20121 and their Figures 3 and 4, for the 
-^ns - -Rns confidence regions). The same data sets are 
used in the present work. This source is moderately 
pilcd-up (~ 4% pile-up fraction) and ne cessitates th e in- 
clusion of a pile-up model component (jDavisl I2001L see 
§ 13.41 for details). All uncertainties for values obtained 
from a X-ray spectral analysis with XSPEC are quoted 
at the 90% confidence level, unless noted otherwise. 

The qLMXB X7 in 47 Tuc has also been observed with 
the high S/N that could provide constraints on the dEoS. 
However, it suffers from a significant amount of pile-up 
(pile-up fraction ~ 10 — 15%). While the effects of pile-up 
can be estimated and corrected for by the inclusion of a 
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pile- up model (|Heinke et al.ll2003l l2006f ). the uncertain- 
ties involved with such a large amount of pile-up are not 
quantified in this model 2 . It was chosen not to include 
this target in the present analysis, in an effort to limit 
the sources of uncertainties that are not quantified (see 

§E3J. 

The list of targets and their usable observations with 
XMM and Chandra is presented in Table [1] along with 
the usable exposure time and other relevant parameters 
of the observations. 

3.2. Data Processing 

The processing of raw data sets is performed according 
to the standard reduction procedures, described briefly 
below. 

3.2.1. Chandra X-ray Observatory Data Sets 

The reduction and analysis of Chandra data sets 
(ACIS-I or ACIS-S) is done using CIAD v4 . 4. The level- 
1 event files were first reprocessed using the public script 
chandrajrcpro which performs the steps recommended 
in the data preparation analysis thread 3 (charge trans- 
fer inefficiency corrections, destreaking, bad pixel re- 
moval, etc, if needed) making use of the latest effective 
area maps, qua ntum efficiency maps and gain maps of 
CALDB v4.4.8 (jGraessle et al.ll2007t) . The newly cre- 
ated level- 2 event files are then systematically checked 
for background flares. Such flares were only found in the 
middle and at the end of an observation of M28 (ObsID 
2683), for a total of 3ks. These two flares caused an 
increase by a factor of 2.4, at most, of the background 
count level. Given the extraction radius chosen here (see 
below), this period of high background contaminates the 
source region with < 1 count. Therefore, the entire ex- 
posure of the ObsID is included in the present analysis. 

This recent version of CALDB contains an improved 
model for the molecular contamination affecting ACIS 
observations. The improvement is noticeably significant 
in the 0.3-0.5 keV range, where the RMS residuals are 
now limited to 10% in that range (Figure 15 of docu- 
ment "Update to ACIS Contamination Model", Jan 8, 
2010, ACIS Calibration Memo 4 ). The 0.3-10.0 keV en- 
ergy range is used for the present analysis. However, 
10%-systematic uncertainties in the 0.3-0.5 keV range 
are included to account for the variations in the resid- 
uals of contamination model. For the remaining energy 
range, 0.5-10 keV, the usual 3% systematic uncertainty 
is added. This is done with the heasof t tool grppha, to 
set the systematic uncertainties in each energy channel. 

3.2.2. XMM-Newton Data Sets 

The reduction of XMM data sets is completed us- 
ing the XMM Science Analysis System vlO.0.0 with 
standard procedures. The command epchain performs 
the preliminary data reduction and creates the event files 
for the pn camera in the 0.4-10.0 keV energy range, with 
3% systematic uncertainties included. The data sets are 
checked for flares and time intervals with large back- 
ground flares are removed. The total usable time (after 



flare removal) for each observation is listed in Table [TJ 
MOS1 and MOS2 data are not used in the present work 
to minimize the effects of cross-calibration uncertainties 
between detectors. 

3.3. Count Extraction 
3.3.1. On Chandra X-ray Observatory Data Sets 

The count extraction of the source and background 
spectra is performed with the task specextract, as well 
as the calculation of the response matrices and ancillary 
response files (RMFs and ARFs). The centroid position 
is chosen using the reported source positions from previ- 
ous works. The extraction radii are chosen to correspond 
to a 99% EEF (Encircled Energy Fraction), and therefore 
depend on the off-axis angle of the targets. For on-axis 
qLMXBs (M13, M28, NGC 6397 and NGC 6304), counts 
within a 3.4" radius are extracted to create the spectrum. 
This ensures that 99% of the enclosed energy fraction at 
1 keV is included 5 . The qLMXB in uj Cen is at a large off- 
axis angle (~ 4.4') which requires an extraction radius of 
6" to contain 99% of the EEF. This is due to the degra- 
dation of the PSF of the Chandra mirror with increasing 
off-axis angle. Background counts were taken from an 
annulus centered around the qLMXB, with an inner ra- 
dius of 5" (9" for ui Cen) and an outer radius of 50". 
Regions surrounding other point sources detected in the 
qLMXB extraction regions or in the background regions 
are also excluded (radius of 5" or more). For NGC 6304, 
the background region is off-centered with respect to the 
source region, to ensure that the background lies on the 
same CCD chip as the source. 

Figures [2^,-e show the regions used to extract the 
counts and create the spectra of each target. When sev- 
eral observations are available for a target, the largest- 
S/N observation was used to create the figure. 

3.3.2. On XMM-Newton Data Sets 

For the three XMM data sets, the extraction method 
was the same at that described above. Only the source 
extraction radii were different and determined using the 
XMM SAS task eregionanalyse which provides the op- 
timum extraction radius that maximizes the S/N given 
the source position and the surrounding background. 
The optimum radius is 19" for ObsID 01122 of w Cen. 
The encircled energy of the source is therefore 79% at 
1.5 keV. 

For the XMM observations of M13, the close proximity 
of a cataclysmic variable (CV) complicates the task. A 
25" extraction radius and a 12.5" exclusion radius for the 
nearby source are used to create the spectra. It ensures 
that 84% of the total energy from the qLMXB at 1.5 keV 
is encircled 6 . 

Similar to the Chandra data, the background is an an- 
nulus around the source, restricted to remain on the same 
CCD chip as the source, rmfgen and arfgen are then 
used to generate the response matrices files (RMF) and 
the ancillary response file (ARF) for each observations. 

3.3.3. Contamination from Nearby Sources 



2 See "The Chandra ABC Guide to Pile-Up" available at 
http: //cxc .harvard.edu/ciao/download/doc/pileup_abc.pdf 5 Chandra Observatory Proposer Guide vl2.0, fig. 6.7, December 

d http://cxc.harvard.edu/ciao/threads/data.html 2009 

4 available at http : \/\/cxc .harvard. edu\/cal\/memos\/contam_memo .pcff From XMM Users Handbook, fig. 3.7, July 2010 
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Fig. 2. — This figure shows the extraction regions for the five qLMXBs. (a) For M28, ObsID 9132, the three nearby sources are excluded, 
limiting contamination to < 1% within the extraction region (see Table[2]l. (b) For NGC 6397, ObsID 7460, no counts from nearby sources 
fall within the extraction region, (c) For NGC 6304, ObsID 11074, the nearby sources are not contaminating the extraction region, (d) 
For M13, Chandra data are on the left, ObsID 7290, and XMM data are on the right, ObsID 0085280301. The nearby CV is excluded. 
(e) For uCen, Chandra data are on the left, ObsID 1519, and XMM data are on the right, ObsID 0112220101. There is no contamination 
from nearby sources. 



As mentioned above, some qLMXBs lie in close prox- 
imity of other contaminating sources. The most evident 
case is that of M13, observed with XMM (Figured;), in 
which part of the counts from a nearby CV still overlaps 
with the qLMXB extraction region even after excluding 
10" around the CV. In M28, three sources in the proxim- 
ity of the qLMXB require parts of the extraction region 
to be excluded (Figure [2k) , with a minor contamination 
from the nearby sources. For the observations of M28 and 
M13, the fraction of contaminating counts present within 
the extraction region of the qLMXBs is estimated using a 
Monte Carlo sampler which draws counts from the radial 
distribution of encircled energy of ACIS or pn at 1 keV 
(see footnotes in previous subsections). Table [2] lists the 
amount of contamination for each observation of M28 
and for the XMM -■pa. observations of M13. The contam- 
ination over the M28 qLMXB region can be neglected, 
since, in the worst case (ObsID 9132) only 47 counts out 
of 6250 are contamination from nearby sources. The CV 



close to the qLMXB in M13 causes a contamination of 
6% and 9% of the counts in each of the XMM spectra, 
or 4% of the total counts available for the qLMXB in 
M13 (XMM and Chandra spectra combined). Overall, 
contamination for nearby sources represents 0.5% of the 
total count number, all sources combined. This contam- 
ination can be safely neglected, since it will not signifi- 
cantly affect the radius measurement. 

3.4. Pile-Up 

Observations of bright X-ray sources may be subject 
to an instrumental effect known as pile-up. When two or 
more photons strike a pixel on an X-ray detector within 
a single time frame (3.24 sec for Chandra- ACIS observa- 
tions and 73.4ms for XMM-pn observations), the pile-up 
effect causes degradation of the PSF and, more impor- 
tantly for the analysis presented here, a degradation of 
the spectral response. Specifically, the recorded energy 
of the event will be the sum of the two (or more) piled- 
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up photon energies. In addition, grade migration (also 
called photon pattern distortion for XMM) also occurs. 
Although a pile-up m odel exists in XSPEC to take into 
account these effects (jDavisI l2001h . it is chosen here to 
restrict the analysis to mildly piled-up observations. 

None the XMM observations of qLMXBs are piled- 
up, given the short duration of a single time frame on 
XMM-pn (73.4ms). Quantitatively, the count rates of 
the qLMXBs in w Cen and M13 (2.6 x 1(T 2 and 2.7x 1CT 2 
counts per seconds, respectively) correspond to ~ 10 ~ 3 
counts per frame. At those rates, the pile-up is negligible. 
The frame time of Chandra-AGIS in full-frame mode, 
however, is significantly longer than that of XMM-pn 
(compensated by the smaller effective area). The Chan- 
dra observations of the qLMXB in M28 are moderately 
piled-up because of a count rate of ~ 0.043 counts per 
seconds (~ 0.14 counts per frame) which corresponds to 
a pile- up fraction of ~ 5% 7 . Such amount of pile- up 
cannot be neglected and is taken into account using the 
pileup model in XSPEC ([Davisl 120011 ). Other Chandra 
observations of the qLMXBs studied in this work have 
smaller count rates which do not necessitate a pile-up 
correction. 

As stated before, the correction of pile-up fractions 
~ 10% and above (such as that of the qLMXB in 47Tuc, 
observed in full- frame mode on Chandra / ACIS) comes 
with unquantified uncertainties. This is tested by sim- 
ulating piled-up spectra of 47Tuc (10-15% pileup frac- 
tion) and M28 (5% pileup fraction) with their respective 
nsatmos parameters, and then by fitting the spectra with 
the nsatmos model without the pile-up component. The 
best-fit radii of each spectrum is affected by systematics: 
- 10% for M28 and - 50% for 47Tuc. The systematic 
error involved with the pile-up of M28 is smaller than the 
measurement error of i?NS m the present analysis, while 
for 47Tuc, the systematic bias caused by pileup is sub- 
stantially larger than the i?NS measurement uncertainty. 
Therefore, the qLMXB in 47Tuc is not used in this anal- 
ysis to avoid introducing a systematic uncertainty which 
may be comparable in size to our total statistical uncer- 
tainty 8 . 

3.5. Spectral Analysis 

The spectral analysis is composed of two parts. In 
the first one, the five targets are analyzed individually 
and the results are compared to previously published re- 
sults. The second part of the analysis in the present 
work pertains to the simultaneous fitting of the targets, 
with a i?NS value common to all five qLMXBs. Prior to 
the discussion of these two parts, the analysis techniques 
common to the two analyses are described. 

3.5.1. Counts Binning, Data Groups and Model Used 

Once the spectra and the respective response files of 
each observation are extracted, the energy channels are 
grouped with a minimum of 20 counts per bin to ensure 
that the Gaussian approximation is valid in each bin. 

7 Chandra Observatory Proposer Guide vl2.0, fig. 6.18, Decem- 
ber 2009 

8 Note that we find in Section l5T2~2l that inclusion of 47Tuc X7 
in this analysis does not significantly affect our best-fit -Rns value; 
nonetheless, we do not include this data, since we cannot estimate 
the effect of including it on our error region. 



TABLE 2 

Count Contamination from Nearby Sources 



Target 


ObsID 


Number of contaminating counts 


Total 




(detector) 


#1 


# 2 


# 3 


sum 


(%) 


M28 


2693 (ACIS) 


0.16 


0.11 


0.16 


0.43 


0.08% 




2684 (ACIS) 





0.15 


0.32 


0.47 


0.08% 




2685 (ACIS) 


0.05 


0.19 


0.17 


0.37 


0.06% 




9132 (ACIS) 


1.7 


38.9 


6.3 


46.9 


0.75% 




9133 (ACIS) 


1.0 


11.7 


2.7 


15.4 


0.65% 


M13 


0085280301 (pn) 


27 






27 


8.7% 




0085280801 (pn) 


34 






34 


6.2% 


Note. 


— The columns 


1", 


«# 2" 


and 3" 


indicate the 


absolute 



numbers of counts from sources falling within the qLMXB extraction region, 
and "sum" is the simple sum of contaminating counts. For M13, there is only 
one nearby source, not fully resolved with XMM. The last column provides 
the amount of contamination as a percent of the total number of counts within 
the qLMXB extraction region. 

For observations with a large number of counts (>2000 
counts, in ObsID 9132 and 9133 of M28 and ObsID 7460 
of NGC 6397), the binning is performed with a minimum 
of 40 counts per bin. In all cases, when the last bin (at 
high energy, up to 10 keV) contains less than 20 counts, 
the events are merged into the previous bin. 

The spectral fitting is performed using the "data 
group" feature of XSPEC. The spectra of each target 
are grouped together, and each group (corresponding 
to each qLMXB) is assigned the same set of parame- 
ters. The spectral model used is the nsatmos model 
( Hcin ke et al.l 12006) . together with Galactic absorption 
taken into account with the multiplicative model wabs. 
The amounts of Galactic absorption, parameterized by 
Nh {Nh,22 in units of 10 22 atoms cm~ 2 , hereafter), are 
fitted during the spectral analysis, and compared to those 
obtained from a HI Galac tic survey from NRAO data 9 
([Dickey fc Lockmanlll990l ). The Nu values used in the 
present analysis are shown in Table [31 The results ob- 
tained with nsatmos are also com pared with the bes t-fit 
results using the model nsagrav dZavlin et al.lH996h for 
completeness . 

As mentioned before, the pileup model ([Davisl 120011 ) 
is necessary for the spectral fitting of M28 spectra. In 
XSPEC, multiple groups cannot be fitted with differ- 
ent models, so a single model is applied to all groups, 
namely, pileup*wabs*nsatmos. For the spectra of the 
qLMXB in M28, the a parameter of the pileup model, 
called "good grade morphing parameter" is left free. The 
frame time parameter is fixed at 3.10 sec. This value 
corresponds to the TIMEDEL parameter of the header 
(3.14104 sec for the observations of M28) where the read- 
out time (41.04 ms) is subtracted. All the other param- 
eters of the pileup model are held fixed at their default 
values, as recommended in the document "The Chan- 
dra ABC Guide to Pile-Up v.2.2" 10 . Since the targets in 
M13, uj Cen, NGC 6397 and NGC 6304 do not require to 
account for pile-up, the time frame for these four groups 
is set to a value small enough so that the pileup model 
has essentially no effect and the a parameters of the four 
non piled-up sources are kept fixed at the default value, 

9 obtained from the HEASARC Nff tools available a t 
http : \/\/heasarc . nasa.gov\/cgi-bin\/Tools\/w3nh\/w3nh.pl 

10 from the Chandra X-ray Sci- 
ence Center (June 2010), available at 

http : \/\/cxc .harvard. edu\/ciao\/download\/doc\/pileup_abc .pdf 



a = 1. A quick test is performed to demonstrate that 
the pileup model with a small frame time has no effect 
on the spectral fit using the non piled-up spectra of the 
qLMXB in NGC 6304. Specifically, the best-fit nsatmos 
parameters and the x 2 -statistic do not change when the 
pileup model (with a frame time of 0.001 sec) is added, 
as expected. 

3.5.2. Individual Targets 

Prior to the spectral analysis, it is important to verify 
that the qLMXBs do not present signs of spectral vari- 
ability. This is done by considering each target individ- 
ually (without the other targets) and by demonstrating 
that the nsatmos spectral model fits adequately, with the 
same parameters, all the observations of the given tar- 
get. More precisely, for each target, all the parameters 
are tied together and we verify that the fit is statistically 
acceptable. 

In Section I4.1[ the results of the spectral fits of in- 
dividual targets are presented. In order to provide the 
full correlation matrix for all the parameters in the fit, a 
Markov-Chain Monte-Carlo (MCMC) simulation is im- 
plemented (described in § [376]) and the resulting posterior 
distributions are used as the best-fit confidence intervals, 
including the Mns _ -Rns confidence region. These spec- 
tral fitting simulations are performed with the Galactic 
absorption parameters Nh left free. This allows us to ob- 
tain best-fit X-ray-measured values of the absorption in 
the direction of each of the targeted GCs. These best-fit 
values are compared to Hl-dcduced values (from neutral 
hydrogen surveys), and are also used for the remainder 
of the work, when Nh is kept fixed. Finally, the results 
of the individual spectral fits are compared to previously 
published results. 

3.5.3. Simultaneous Spectral Analysis of the Five Targets 

The main goal of this paper is the simultaneous spec- 
tral fit of five qLMXBs assuming an i?NS common to all 
qLMXBs. Therefore, i?NS is a free parameter constrained 
to be the same for all data sets, while each NS targeted 
has its own free Mns and kT c R parameters. In addition, 
the spectra of M28 require an extra free parameter a for 
the modeling of pile-up. This leads to a total of 12 free 
parameters. 

In an effort to include all possible uncertainties in the 
production of the Mns - -Rns confidence regions, Gaussian 
Bayesian priors for the source distances parameters are 
included, instead of keeping the parameter values fixed 11 . 
Since, additional systematic uncertainties can arise when 
keeping the Nh parameters fixed, this assumption is also 
relaxed in the spectral analysis. Finally, the spectra of 
some qLMXBs display excess flux above 2 keV, which is 
not due to the H-atmosphere thermal emission. This is 
accounted for by adding a PL spectral component to the 
model, where the photon index is fixed at T = 1.0 but 
the PL normalizations are free to vary. Such PL index is 
t he hardest observed for a LMXB in quiescence (Cen X- 
4 ICackett et al.l I2010D . The spectra resulting from this 
analysis (with all five qLMXBs, and with all assumptions 
relaxed) are shown in Figured! 

Relaxing the assumptions mentioned above adds 15 
free parameters, for a total of 27, which increases the 

11 For a review on Bayesian analysis, see I Gregory! (120051 ') 



complexity of the x 2 -space. Because of that, in XSPEC, 
the estimation of the confidence region for each param- 
eter proves difficult. The command steppar iteratively 
calculates the \ 2 -value for fixed values of a parameter in 
the range provided by the user. However, this grid-search 
procedure is highly dependent on the starting point of the 
parameter of interest and on the number of steps. Such 
a problem is particularly evident in the case of highly co- 
variant sets of parameters. This can result in ID or 2D 
Ax 2 contours that are not reliable to estimate the uncer- 
tainties. A solution to this issue consists of using the pos- 
terior distributions from MCMC simulations (described 
in § I3.6|) . With those, one can quantify the uncertainties 
of each parameters. The need to include Bayesian priors 
also brings forward the use of MCMC simulations. 

The simultaneous spectral fitting of all five targets us- 
ing MCMC simulations was performed in seven sepa- 
rate runs, during which the assumptions on the spectral 
model are progressively relaxed. The characteristics of 
each run are described in Section 14.21 Another run is 
also performed with the nsagrav model for purposes of 
comparison with the nsatmos model. The following sub- 
section describes the MCMC analysis performed. 

3.6. Markov- Chain Monte Carlo Analysis 

As described above, the main advantage of using an 
MCMC simulation resides in a complete understanding 
of the posterior probability density functions of each pa- 
rameter. It also allows one to marginalize over the so- 
called nuisance parameters, i.e., those that are an impor- 
tant part of the modeling but which are of little physical 
interest to the problem at hand. 

Due to the skewed nature of the parameter spaces, 
leading to curved Mns _ -Rns contours, an MCMC al- 
gorithm different from the typically used Metropolis- 
Hasting algorithm (MH) is chosen. Indeed, we find that 
the MH algorithm is not efficient at exploring skewed pa- 
rameter spaces. The next few paragraphs are dedicated 
to a brief description of the Stretch- Move algorithm used. 

Th e Stretch-Move algorithm (|Goodman fc Weard 
I2010T) is particularly useful for elongated and curved dis- 
tributions (e.g., Mns~-Rns wit h nsatmos), as de mon- 
strated in previous works fe.g.. IWang et all 120111 ); our 
implementation generally follows that work. Other anal- 
yses have used the Stretch - Move algorithm (jBovv et al.1 
I20T1 lOlofsson et all I20TI using ano ther implementa- 
tion, by iForeman-Mackev et aill2012f) . The algorithm 
consists of running several simultaneous chains, also 
called walkers, where the next iteration for each chain 
is chosen along the line connecting the current point of 
the chain and the current point of another randomly se- 
lected chain. The amount of "stretching" is defined by 
a random number z from an affine-invariant distribution 
(|Goodman fc Weardl20Tol) : 



1 if ze [i ,ol 
otherwise 



(2) 



In this work, each individual chain starts from a ran- 
domly selected point in the parameter space, within the 
typical hard limits defined for nsatmos in XSPEC. 

The parameter a is used as a scaling factor and is ad- 
justed to improve performance. A larger value of a in- 
creases the "stretching" of the chains, i.e., the algorithm 
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Fig. 3. — This figure shows the spectra resulting from Run #7, obtained with the model wabs* (nsatmos+pow) , for which the best-fit 
statistic is xg/dof (prob.) = 0.98/628 (0.64). 



will better explore the elongated parts of the parameter 
space, but it decreases the likelihood that the next step is 
accepted. A smaller value of a only produces small excur- 
sions from the previous value but increases the likelihood 
that the next step is accepted. Efficiency is optimized at 
an intermediate value of a. The Stretch-Move algorithm 
can be fine-tuned with only two parameters: a and the 
number of simultaneous chain. By comparison, the MH 
algorithm requires N(N + l)/2 tuning variables, where 
N is the number of free parameters. 

The validity of this MCMC algorithm is assessed 
by performed a test run with a single source (U24 in 
NGC 6397, with fixed distance) and comparing the re- 
sulting Mns - -Rns contours with those obtained from a 
simple grid-search method (steppar in XSPEC). Specif- 
ically, the obtained Mns - -Rns contours as well as other 
posterior distributions match those obtained from a 
steppar grid-search in XSPEC. The addition of Gaus- 
sian Bayesian priors on the distance is also tested with 
U24, which results in Mns~-Rns contours broadened in 
the i?oo direction. This is because the normalization of 
the thermal spectrum is approximately cx (Roo/d) 2 . 

For the Stretch-Move algorithm, the minimum num- 
ber of simultaneous chains is equal to N + 1, where N is 
the number of free parameters. However, increasing the 
number of simultaneous chains ensures a more complete 
coverage of the parameter space, when comparing the re- 
sults of the Stretch-Move algorithm to contours obtained 
with steppar. In addition, it reduces the chances of hav- 
ing the N + 1 walkers collapsing to a N - 1 dimensional 



space, i.e., one of the parameters has the same value 
within all the chains causing all following steps to evolve 
in the same plane. However, increasing the number of 
walkers also increases the convergence time. 

The resulting posterior distributions are then 
marginalized over nuisance parameters. While necessary 
for the spectral fitting, these parameters do not provide 
physical information (e.g., a, the pile-up parameter). 
The results are presented in Section @J where the values 
quoted correspond to the the median value (i.e., 50% 
quantile) of each parameter. For highly skewed param- 
eter posterior distributions, the median value can be 
very different from the most probable value. The error 
regions (90%-confidence) of each quoted parameters 
are also calculated with quantiles. In some cases, the 
normalized probability of a parameter posterior distri- 
bution does not converge to zero within the parameter's 
hard limits in XSPEC. This is indicated by a 'p' in the 
tables listing the parameters. 

3.7. Distances to the Globular Clusters and their 
Uncertainties 

While most GCs have distances estim ated from pho- 
tome tr y - using RR Lyrae variable stars (jMarconi et al.l 
120031: IBono et alj 120071). horizontal branch stars 
(jValenti et alj|2007t IGratton et al.ll2010D . or the carbon- 
oxyge n white-dwarf (CO- WD) sequence ((Hansen et al.1 
l2007f ) - these methods suffer from systematic uncertain- 
ties that are difficult to quantify. In fact, many recent 
photometric studies of GCs do not quote the amount of 
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TABLE 3 

Globular Cluster Relevant Parameters 



Name cfcc (kpc) Method Nh.22 (X-ray) Nh 22 (HI) Reference 

M28 5.5±0.3 Horizontal Branch fitting 0.256^^^ 0.24 Testa et al. (20011 

NGC 6397 2.02±0.18 Dynamical 0.096±j};[J£[ 0.14 Rees (19961 

M13 6.5±0.6 Dynamical 0.008to.oo7 0.011 Rees (19961 

uCen 4.8±0.3 Dynamical 0.182±°'°g 0.09 van de Ven et al. (20061 

NGC 6304 6.22±0.26 Horizontal Branch fitting 0.346±o HI 0.266 Recio-Blanco et al. (20051 



Note. — The selection of the distance values is described in Scction l3.7l and the quoted uncertainties are la. The 
Nh values are given in units of 10 22 atoms cm -2 , with 90%-confidcncc uncertai nties from X-ray spectral fitting. The 
Nh (HI) column corresponds to value in the direction of GCs, in the HI survey of (Dickey & Lockman 1990). The X-ray 
values arc deduced from the best-fit N h obtained from X-ray spectral fitting of each target in this work. Only the N h 
values for NGC 6397 and uj Cen are not consistent with the HI values (sec § 14. II for details). Nh values deduced from 
the present X-ray spectral analysis are used in the present work. 



uncertainty in the measured distance (iRqsenberg et al.1 
2000; IBica et al.l [200(1 IGratton et al.ll2010D . 

While some references discuss systematic uncer- 
tainties related to th e correction of extinction (e.g. 
iMcDonald et al.l [2009L for wCen), other sources of sys- 
tematic errors can affect the results, including errors re- 
lated to the m etallicity of c luster members (see dispersion 
in Figure 1 of Harris 20j3), to a possible differential red- 
dening in the d irection of GCs (as observed for oj Cen, 
lLaw et al.ll200l . to variations in the modeling of extinc- 
tion w ith R(V) ~ 3.1 — 3.6 ( adding ~ 10% of uncer- 
tainty, iGrebel fc Roberts) 1 1995( 1. or to the stellar evolu- 
tion/atmosphere models used. As an example for the lat- 
ter, distance determination methods to NGC 6397 using 
CO- WD may be affe cted by uncertainties in the evolu- 
tiona ry code models (|Hansen et al J 120071 : iStrickler et al.l 
I2009T ). which are not easily quantifiable. 

Therefore, whenever possible, dynamical distance mea- 
surements are used - distances estimated from proper 
motion and radial velocities of cluster members. These 
purely geometrical methods produce well-understood un- 
certainties, although they are at the moment larger than 
reported uncertainties from photometric methods. This 
is consistent with the goal of this paper which is to es- 
timate i?NS and its uncertainties, minimizing systematic 
uncertainties. The upcoming mission GAIA from the 
European Space Agency, scheduled for 2013, is expected 
to produce GC distance measurements, to an accuracy of 
few p ercent, by determining the parallax of cluster m em- 
bers (jBaumgardt fc KroupalfeOOa lBaumgardtll2008f ). 

The adopted distance values are discussed below and 
are summarized in Table [3l In the following list, un- 
certainties are quoted at the la level (for GC distances, 
distance modulii, etc.) 

• The GC M28 does not have a dynamical distance 
measurement, but its distance has been estimated 
in dif ferent works: 5.1 ±0.5kpc dRees fc Cudworthl 
[l99lh . 4.8-5.0 kpc (IDavidee et alJ 11996ft and 
5.5 kpc (|Harrisl[T99l iTesta et al.ll2001l )~ all using 
photometric methods. For the most recent re- 
sult, uncertainties can be estimated from the uncer- 
tainties in the horizontal branch (HB) magnitude. 
Specifically, the uncertainty in Vhb = 15.55 ± 0.1, 
translates into the u ncertainty in the di stance: 
d M2S = 5.5 ± 0.3 kpc (jServillat et al.|[2012T ). This 
measured value and its uncertainties were used 
here. 



• The distance to NGC 6397 has been reported 
from a dyna mical study to be dNGC 6397 = 
2.02±0.18kpc (|Reesl 119961) . More recent photo- 
metric studies (CO WD sequence) h ave been per- 
form ed, with d = 2.54±0. 07kpc (|Hansen et al.l 
120071 ). or d = 2.34±0.13kpc (|Strickler et al.1 12009). 
but since those results are model-dependent, they 
are not used in an effort to minimize unquantified 
systematics. When the present analysis was at an 
advanced near-completion stage, recent results re- 
porting a dynamical measurement of the distance 
came to our attention: ^ncc 6397 = 2.2^7 kpc 
(|Hevl et al.l l2012f) . consistent with g?ngc 6397 = 
2.02±0.18 kpc, the value used in the present work. 

• For M13, the dynamical di stance has been mea- 
sured: d M13 = 6.5±0.6kpc (|Reeslll996h . No other 
paper in the literature reports a distance mea- 
surement with quantified uncertainty. This dy- 
namical measurement is consistent with th e value 
divri3 = 7.1 kpc obtained for m photometry ([Harris! 
1996: iSandquist et alJ[201f3 ). While uncertainties 
could be estimated for this measurement like it was 
done for M28, the dynamical measurement is pre- 
ferred to limit the effect of systematic uncertainties, 
as explained above. 

• u> Cen's distance was measured in a dynamical 
study , rf w Cen = 4.8±0.3kpc (jvan de Ven et al.l 
2006), and no other reference provides a distance 
with its measurement uncertainty. This mea- 
surement is consi stent with o ther estimates (e.g. 
d wCen = 5.2 kpc. [Harriill996l update 2010). 

• The GC NGC 6304 lacks a dynamical distance 
measureme nt. However, results fr om a previ- 
ous work (jRecio-Blanco et all 120051. using pho- 
tometric data from iPiotto et alJ I2002T ) are avail- 
able. In that work, the distance modulus in the 
F555W filter (Hubble Space Telescope filter) is 
(to - M) F555W = 15.58±0.09. The reddening in 
this band for NGC 6304 was not provided in the 
published work, bu t the value E(B — V) = 0.52 
(jPiotto et all [2002) can be used instead. This is 
acceptable because the averag e difference between 
E (HS T) and E (B - V) in the iRecio-Blanco et all 
(2005) catalogue is AE = 0.005, which has a 
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TABLE 4 

Spectral Fit Results of Individual Sources 



Target 


fcT cff 
(keV) 


-Rns 


M NS 


Roc N H ,22 


xl/d.o.l. (prob.) 




(km) 


(M ) 


(km) 





M28 120^*2 iO^^j.g 1 - 25 -o.li P 13 -°li'.9 °- 252 la024 094 / 269 (°- 76 ) 

NGC 6397 76±" 6.6l^ 0.84±g |g p 8.41J - ? 0.0961q 1.06 / 223 (0.25) 

M13 83±H W.1+1% 1-271°;™, 12.8+H O.OQSl ,^ 0.94 / 63 (0.62) 

wCen 64l? 7 20.1+^ p 1.781 J;^ 23.61?;? 0.182l°'°« 0.83 / 50 (0.80) 

NGC 6304 1071?? S^l^L 1.16±g;gg p 12-gl|g °- 346 1q.093 1 - 07 / 29 (°- 36 ) 

M28 119lg a 10.6^2.6 1 - 17 -o.56p 12 - 9 lo'.9 (°- 252 ) 0.94 / 270 (0.77) 

NGC 6397 76lg 5 6.6l?;L 0.84l° j| p 8A±°[ ;f (0.096) 1.06 / 224 (0.26) 

M13 86l?Q 9- 2 t\ 7 Zp i-iSlo.^p 11.6±i | (0.008) 0.93 / 64 (0.63) 

wCen 64l^ 3 19.6l|;| 1.84l?;?^ 23.2±|;| (0.182) 0.82 / 51 (0.82) 



_ 5 ^-"-3.8 ^'^-l.lOp * J -*-3.3 

NGC 6304 106l?3 9A -2A P ^lo'Jip 1L8 l2.o (0.346) 1.05 / 30 (0.39) 



Note. — The targets were fit individually with fixed distances and free/fixed Nh values. For 
M28, the pileup model is included (see § 13.51 for details), and a value a — 0.45^q ^| is obtained. 
The posterior distribution of Hoc was obtained by calculating the value of i?^ from -Rns and A/ns 
at each accepted MCMC iteration. Quoted uncertainties arc 90% confidence, "p" indicates that 
the posterior distribution did not converge to zero probability within the hard limit of the model. 



negligible effect on the absolute distance modu- 
lus. Therefore, (m - M) = 13.97±0.09, assum- 
ing Ay = 3.1 x E(B-V), give d NGC 6 304 = 
6.22±0.26kpc. 

Overall, the distances to the targeted GCs have un- 
certainties of ~ 9% or less, keeping in mind that the 
distances determined with photometric methods possi- 
bly have systematically underestimated uncertainties. 

4. RESULTS 

In this section, the results of the spectral analyses of 
each target individually, with their measurements, 
are first presented. These include comparisons with pre- 
viously published results. In particular, some issues re- 
garding the reported spectral analyses for the qLMXBs 
in oj Cen and M13 are raised. Following this, the results 
of the i?NS measurement from the simultaneous fit are 
detailed. 

4.1. Roo Measurements of Individual qLMXBs 

The analysis of the targeted qLMXBs is performed 
with the spectral model detailed above (§ I3.5[) . For each 
target, analyzed individually, the fits are statistically ac- 
ceptable (i.e., with a null hypothesis probability larger 
than 1%), which demonstrates that, within the statis- 
tics of the observations, the sources did not experience 
any significant spectral variability over the time scale be- 
tween the observations. The resulting values and 90% 
confidence uncertainties, along with the x 2 -statistic ob- 
tained, are provided in Table [4] The spectral results 
obtained with Nh fixed at the X-ray-deduced values, in- 
stead of the usual HI survey values, are also provided. 
Discrepancies between the X-ray-deduced and HI survey 
values of Nh, if any, are discussed for each individual 
target. 

Table [4] also shows the best-fit Roo values, calculated 
using the equation: 

(2GAI \ ' " 



from each accepted points of the MCMC runs. Uncer- 
tainties in Roo are then obtained from the calculated 
posterior distributions of Roo resulting for the MCMC 
runs. The use of MCMC simulations has the advantage 
of avoiding geometrical construction to calculate the un- 
certainties icJ_J?2^Jroin the M ns-Rns contours as per- 
formed in lGuillot et all (|2011aT l. 

In the following subsections, the previously published 
results are compared to those obtained here. To do so, 
the Roo measurements are renormalized to the distance 
used in the present analysis. 

4.1.1. Comparison with Published Results - M28 

Using the 2002 Chandra data, the reported 
Roo value of the qLM XB in M28 was Rqq = 
14.5l|| km (£>/5.5kpc) (IBecker et all I2003I ). An 
additional 200 ks of observations obtained with Chan- 
dra-ACIS in 2008 was used to produce a refined radius 
measurement: i?NS = 9±3km and Mns = 1.4l °; 9 Mq , 
with an H- atmosphere model (jServillat et al.1 I2012T ). 
corresponding to Roo — 12.2l 1 ; 4 km for dM28 = 5.5 kpc, 
consistent with the discovery work (jBecker et a,l.ll2003h . 

All the NS parameters resulting from the present 
analysis (Table H Roo = 13.0l?;^km (£>/5.5kpc), 

for Nh,22 — 0-2521qq24) are a l so consistent with the 
previously published results. In addition, the previ- 
ous work also perfor med a careful variability analysis 
(jServillat et al.ll2012| ). confirming our findings that the 
qLMXB in M28 is not variable. 

The best-fit Nh found here is consi stent with the value 
from an HI survey: Nh,22 = 0.24 (Dic key fc Lockmanl 
1990), but the X-ray- measured Nh value is preferred in 
the rest of the present work, for the MCMC runs with 
fixed Nh- 

4.1.2. Comparison with Published Results - NGC 6397 

The data sets used in this wor k are the same as th e 
ones used in the previous work (Guill oFet al.l l2011a|) . 
There are however minor differences in the data reduc- 
tion, namely, the extraction radius used (99% EEF in this 
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Fig. 4. — This figure shows the one- and two- dimensional marginalized posterior distributions for the NS properties (radius, temperature 
and mass) obtained from the MCMC run of the qLMXB in M28. The ID posterior probability distributions are normalized to a total 
probability of unity. The top-right plot shows the ID posterior distribution of Roo values. The 68%, 90% and 99% confidence intervals or 
regions arc shown with solid, dashed, and dotted lines, respectively. In each of the ID distributions, the median value is shown as a red 
line. Note that the 99% region is not always visible in the ID distributions. The physical radius of the NS in M28 is -Rims = 10.6^2 6 
This corresponds to a projected radius of Roo = 12.9^g'g km, for Njj,22 = 0.252. The double-peaked ID distribution of Mns is due to 
the strongly curved nature of the Mns — -^NS an d A^NS - kT e g 2D distributions, i.e., the strong cor relation betw een these parameters. This 
figure and the following Figures IBT I16I were created with the Mathematica package LevelSchcmes ( Caprio 2005|) . 
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Fig. 5. — This figure is similar to the precedent (Fig. gj, but for the qLMXB in NGC 6397. The physical radius of the NS is R NS = 
6-6l5'Ip km which corresponds to = 8-4+g"| km, for N Ha 2 = 0.096. 



work compared to 98% EEF at 1 keV previously), the cal- 
ibration files used (latest version of CALDB v4.4.8), the 
distance used for the spectral fit, and the energy range 
(0.5-8 keV in lGuillot et alll2011a|) . 

After re-normalizing to the distance used 
in the present work, the previous R^ result, 
= 9.6±£|km (L>/2.02kpc), is consistent 
with the one obtained from the MCMC run: 
i?oc = 8.4t^km ( J D/2.02kpc), for N H ,22 = 0.096±g;g^. 
This best-fit value of Nh is howeve r, inconsistent with 
the fixed HI value (Nh = 0.14, Dickcv fc Lockmanl 
1990]) used in the previous work dGuillot et al.ll2011al) 12 . 

12 The X-ray deduced value of Nu found here is nonetheless 



This puts into question the i?NS measurement and 
-Mns - -Rns contour s previously publis hed with the 
value iV ff ,22=0.14 (IGuillot et all 1201 lafl . When fixing 
Nh, 22=0-14: in the present work, the resulting i?^ value 

is Roo — 11.9 loskm (J/2 .02 kpc), marginally consis- 
tent with the (jGuillot et alll2011al ) result. Nonetheless, 
one notices that the different value of Nh causes a 
significantly different resulting i?NS value. Basically, 

con sistent with the Nh v alue from a different survey of Galactic 
HI UKalberla et al.ll2005t) . ##,22=0.11, and with the N H value 
calcu lated from the reddening in the direction of NGC 6397 l|Harrisl 
119961) with a linear relati on between Njj and the extinction Ay 
UPredehl fc Schmittiri995D . 
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Fig. 6.— This figure is similar to Fig. g] but for the qLMXB in M13. The physical radius of the NS is R NS = 9- 2 ±2.3 P km which 
corresponds to R^ = ll.e+^J km, for N Ha 2 = 0.008. 



increasing the assumed value of Nh for a given target 
leads to a larger R^. This is further discussed in 
Section [5j In the rest of the present work, the best-fit 
X-ray deduced Nh value Nh,22 = 0.096 is used. 

4.1.3. Comparison with Published Results - oj Cen 

The original Roo measurement from the Chandra dis- 
covery observations wjs_R 2SL ^14 ; 3±2JJtin (D/5.0kpc) 
for N H ,22 = 0.09 (jRutledge et al.l I2001 . or R^ = 
13.7±2.0km (D/4.8 kpc). Another work measured 
Roo = 13.6±0.3km (D/5.3kpc) with N H22 = 0.09 ± 
0.025, equivalent to R x = 12.3±0.3 km (D/4.8 kpc) 
using the XMM observation of u) Cen (|Gendre et al.l 
l2003a|) . Results from both analyses are consistent with 



the radius measurement performed in this work, with 
the value of Galactic absorption 7X^22 =0.09: Roo — 
11.9t^km (D/4.8 kpc). However, when removing the 
constraint on Njj, the best-fit Roo and Njj become incon- 
sistent with the previously reported values. Specifically, 
Roo = 23.6±£?km (D/4.8 kpc) for N H , 2 2 = 0.182±g;^. 
This value of Nh , not co nsistent with the HI survey value 
(jDickev fc Lockmanlfl99l . was used in the remainder of 
the present work. One can also note that the present re- 
sults (best-fit i?NSj -^NSj kT c g, and Nh) are consisten t 
with those previously published (Webb fc Barretll2007[ ). 

The results presented in Table |4] should be treated 
as more realistic than the initially reported one since 
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Fig. 7. — This figure is similar to Fig. [4] but for the qLMXB in o;Cen. The physical radius of the NS is -Rns 
corresponds to = 23.2^30 km, for N Ha2 = 0.182. 



lg.eig j km which 



they make use of more recent calibrations of XMM and 
Chandra, as well as an improved method. In particu- 
lar, the s mall uncertainties (~ 2%) on i?^ previously 
published (jGendre et alJl200"3al ) are particularly intrigu- 
ing. It has also been shown in another reference that the 
S/N obtained with 50 ks exposure of uj Cen is not suf- 
ficient to constrain th e radius with ~ 2% uncertainty 
()Webb fc Barreti 120071 ). but the cause of this discrep- 
ancy was not discussed. The co nstrained J?nr mea- 
surement with ~ 2% uncertainties (IGendre et al.ll2 003a') 



was not reproduced in the later work (jWebb fc Barreti 
120071 R m = 11.7i^km, using the same XMM data), 
nor in the present work. Using the same model as the 



one initially used (|Gendre et al.l [2003a|) . similar uncer- 
tainties (~ 2%) can only be obtained when keeping the 
NS surface temperature fixed, leaving the normalization 
(i.e., the projected radius i?oo) as the sole free param- 
eter. Specifically, with the same model and analysis 
procedure, the uncertainties on i?^ are aj^ ~ 3% 

with the temperature fixed and becomes aj£ ~ 15% 
when the tempe rature is a free param eter. If this is the 
method used in IGendre et al.l (|2003aD . the uncertainties 
of i?oo only represent the statistical uncertainties and 
are therefore highly underestimated. It is inappropri- 
ate to keep the temperature fixed because there is no 
known prior on the NS surface temperature, and there- 
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Fig. 8.— This figure is similar to Fig. [4] but for the qLMXB in NGC 6304. The physical radius of the NS is i? NS = 9.5jl|;l, km. The 
corresponding projected radius is = ll.sj^'rj km> for Nn,22 — 0.346. 



fore it must remain free during the spectral fitting. In 
addition, the XMM-pn observations suffer from periods 
of high-background activity which need to be removed 
(see Figure [T]). This leads to 34 ks of usable exposure 
time of the 41 ks available. No su ch background flares 
were reported in the o riginal works ( Gen dre et al.l feOOSa; 
iWebb fc Ba7retl[2007t l. 

This note about the amount of uncertainty for w Cen 
is of crucial importance since this source has often been 
cited as the canonical qLMXB, with the best radius 
measurement available, citing the underestimated ~ 2% 
uncertainties on i?NS ■ Deeper exposures of u> Cen are 
needed to provide constraints that will be useful for dEoS 



determination. Moreover, this discussion also points out 
the importance of reporting Mns~Rns contours (instead 
of simple -Rns measurements) for the measurements of 
NS properties using the thermal emission from qLMXBs. 

4.1.4. Comparison with Published Results - M13 

The Roo value of the qLMXB in M1 3 reported in the 
discovery paper (|Gendre et all I2003bfl . fi M = 12.8 ± 
0.4km (D/7.7kpc), corresponds to = 10.8 ± 

0.3km (D/6.5 kpc). This value is consistent with the 
value presented in the present work, given the uncer- 
tainties: = 12.8l|;£km (D/6.5kpc), for N H ^2 = 
0.008±^o7 P - The best fit N H is consistent with HI sur- 
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TABLE 5 

Results from Simultaneous Spectral Fitting, with Fixed N h 

Target a p ii eU p kT cS M m N H ,22 PLNormxlO- 7 

(eV) (M ) (km) keV" 1 s" 1 cm" 2 



Run #1: Fixed N H , Fixed d GC , No PL included, _R NS =7.1_gg km 

xj/dol (prob.) = 97/643 (0-70), 18% accept, rate 

M28 176+it 1.62+^'^ (0.252) ~ 



NGC 6397 71± 7 3 0.69t°;^ 8.4j:°-j (0.096) 

M13 110 ±io 1 - 41 i u ;29 n -°-i.i (0-008) 



uCen !64+i4 2.05t°" 18. 9±l'l (0.182) 

NGC 6304 - 136±^ 1.4ltp.43 P jj^jig (0-346) 

Run #2: Fixed Nh , Gaussian Bayesian priors for d G c , No PL included, -RNS=7.6^Qg 1 

X^/dof (prob.) = 0.98/638 (0.64), 11% accept, rate 

M28 0.44+^;^ 165+^ 1.63+JJ;i| 12.6^; (0.252) 
NGC 6397 71±| 0.73±g |g p 9.0tJ; \ (0.096) 

M13 iOltig 1 -34to.53 H-Oi};? (0.008) 

i^Cen 154^22 2 - 16 ±o.2l 19 - St^l (0-182) 

NGC 6304 - 127l 2 ^ ^gtpj^ H-Oti's (0-346) 

Run #3: Fixed N H , Fixed d GC , PL included, i? NS =7.3±og km 

Xg/dof (prob.) = 0.96/638 (0-78), 15% accept, rate 

M28 0.35iJJ;^ 170^ 1.63^;^| 12. 6+^ (0.252) 5.li^;^ p 
NGC 6397 70tg 0.68±°' 2 g 8.6±og (0.096) 2 - 2 ±i'.3 p 
M13 109+J 2 1.50l°;^f 11-7±J;| (0.008) 2 - 2 ±2.o P 
a; Cen 163t^ 2.14^0-14 19.9±^ (0.182) 2.9±| J p 
NGC 6304 136 ±^ l-51±g || 1P 11.7± 2 ;[j (0.346) i-Tti'jp 



Run #4: Fixed TVjf , Gaussian Bayesian priors for d GG • PL included, -Rns=8.0^ j'JJ ] 

X 2 /dot (prob.) = 97/633 (0-72) 11% accept, rate ' 

M28 0.35l[j;^ 157+^ 1.64+^;J| 12. 8^ (0.252) 5-0 -tip 

NGC 6397 70^3 0.72±°H 9A+ (0.096) 2.2+ 1 ' 3 



M13 99± 22 l-43to6L n - 7 -ii (0.008) 2.2±^° 

1+21 na +6.2i on .+2S / n1 „,i , n+3-8 



Cen 1511'i 2 - 2 8^.^ 20.4^° (0.182) 3.0l 2 . 5p 

NGC 6304 l 2 5tao 1A6 -o:ltp n - 8 -2 2 i (0.346) 



Note. — a p ii cup corresponds to the parameter of the pileup model. "PL Norm." refers to the value of 
the normalization of the power-law component, when used. For each run, the characteristics arc described: 
whether or not the absorption Nh was fixed; whether the GC distances dec were fixed or if a Bayesian 
prior was imposed; whether or not a additional power-law component (PL) was included in the model. For 
each run, the best xt value is provided, as well as the null hypothesis probability. Finally, the acceptance 
rate (not including the burn-in period) is provided. All quoted uncertainties arc 90% confidence. Values 
in parentheses arc kept fixed in the analysis, "p" indicates that the posterior distribution did not converge 
to zero probability within the hard limit of the model. 



vey values (N H = 0.011, iDickev fc Lockmanlll990l) , but 
Nh,22 = 0.008 is used in the remainder of the present 
analysis. 

Once again, the uncertainties reported in the original 
work are small and can only be reproduced when fixing 
the temperature. Similarly to u> Cen, the M13 discovery 
analysis was likely performed keeping the temperature 
frozen to estimate the uncerta inty on the radius, and 
the ±0.3 km uncertainties cited (jGendre et al.ll2003bD are 
only systematic uncertainties. 

In summary, the results found in this work for M13 are 
consistent with the existing ones (|Gendre et al.l I2003bt 
I Webb & Barret] 120071 ). and while our radius measure- 
ment uncertainties are not as constraining as those previ- 
ously reported, they are considered more realistic given 
the S/N available for the observations, and given that 
more recent calibrations have been used. Similarly to 
uj Cen, deeper exposures of M13 would provide the nec- 
essary S/N to constrain the dEoS. 



When the present work was at an advanced stage, re- 
sults of a spect ral analysis of the qLM XB in M13 came to 
our attention (jCatuneanu et al1 l2013). These results are 
consistent with those found in the present work, when 
re-normalized to the distance used here. 

4.1.5. Comparison with published results - NGC 6304 

This analysis presents a new 100 ks observation 
of NGC 6304. The Roo value of the qLMXB, 
Rao = 12.2^ km (D/6.22kpc), for N H . 22 = 
0.3461q'q93), is consistent with that obtained from 
the XMM observation (|Guillot et al.l l2009a|) . R^ = 
12.1±4gkm (£>/6.22kp c), and with t hat fro m a short 
Chandra observation (jGuillot et al.1 12009b). R^ = 
7.8l3;gkm (D/6.22kpc), after re-normalizing the 2009 
measurement to the distance used in the present paper. 
The best-fit X-ray deduced value is also consistent with 
the value used in the original work, obtained from HI 
surveys. Nonetheless, the X-ray measured Nh value is 
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Fig. 9. — This figure shows the marginalized posterior distribution in JW^s - -Rns space for the five qLMXBs, in the first MCMC run, 
where the distance and the hydrogen column density Nu are fixed and where no PL component is added, corresponding to Run #1. The 
68%, 90% and 99%-confidence contours are shown with solid, dashed and dotted lines on the M^s-R^s density plots, respectively. The 
bottom-right graph is the resulting normalized probability distribution of i?NS> common to the five qLMXBs, with the 68%, 90% and 
99%-confidence regions represented by the solid, dashed and dotted vertical lines. The median value is shown by the thin red line. The 
measured radius is i?Ng = 7.1^Qg km (90% confidence). 

used in the remainder of the present work. 

4.2. i?NS; Measurement of the Radius of Neutron Stars 
In this section, the results of the simultaneous spectral 



fits with the parameter posterior distributions obtained 
from the MCMC simulations are presented. The follow- 
ing distinct MCMC runs are performed: 

• Run #1: Model nsatmos with fixed Nh values 
and fixed distances: 12 free parameters, 25 Stretch- 
Move walkers. 

• Run #2: Model nsatmos with fixed Nh values and 



Gaussian Bayesian priors for the distances: 17 free 
parameters, 30 Stretch-Move walkers. 

• Run #3: Model nsatmos with fixed Nh values and 
fixed distances, and an additional PL component: 
17 free parameters, 30 Stretch-Move walkers. 

• Run #4: Model nsatmos with fixed Nh values, 
Gaussian Bayesian priors for the distances, and an 
additional PL component (with fixed index L = 
1.0, but free normalizations): 22 free parameters, 
35 Stretch-Move walkers. 
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Fig. 10. — This figure is similar to the previous one, Fig [9] but for 
the distances to the five qLMXBs (see Table Q5). The resulting radius 

• Run #5: Model nsatmos with free Nh values and 
fixed distances: 17 free parameters, 30 Stretch- 
Move walkers. 

• Run #6: Model nsatmos with free Njj values and 
Gausssian Bayesian priors for the distances: 22 free 
parameters, 35 Stretch-Move walkers. 

• Run ^=7: Model nsatmos with free Nh values, 
Gaussian Bayesian priors for the distances, and an 
additional PL component: 27 free parameters, 40 
Stretch-Move walkers. The spectra resulting from 
this run are shown in Figure [3J 

• Run ^8: Model nsagrav with fixed Nh values 
and fixed distances: 12 free parameters, 25 Stretch- 
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the MCMC Run #2, where Gaussian Bayesian priors were used for 
measurement is i?Ng = 7.6j^'g km. 

Move walkers. This model is used for comparison 
with the nsatmos model. 

All the runs converged to a statistically acceptable point 
in the parameter space, with \t ~ 1 and a null hypoth- 
esis probability > 0.01. In addition, the acceptance rate 
of each run is large enough (> 5%) that the model used 
and the assumptions are adequate for the data. From 
the accepted steps of each run, the marginalized poste- 
rior distributions of all parameters and the median values 
with 90% confidence regions are quoted in the tables. In 
addition, Mns _ -Rns posterior distributions for each of 
the five qLMXBs are obtained. 

Run #1 was performed with the maximum constraints 
imposed on the model. With the assumptions imposed 
on the parameters, this run leads to the most con- 
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Fig. 11. — This figure is similar to Figure [9] corresponding to the results of Run #3. In this run, the distances are fixed (no priors 
included), but a PL component (with fixed index T = 1) is added to the spectral model, leading to -Rns = 7.3jl ' 8 km. 



strained M^s~Rns contours of this work, resulting in 
i?NS = 7.1+^km (90%-confidence), and the minimum 
X^/dof (prob.) = 0.97/643 (0.70). Detailed information 
for other parameters is shown in Table [5] 

The possible effect of auto-correlation between the 
steps of the MCMC simulation is investigated by select- 
ing ev ery other 10 accep ted point, a method called thin- 
ning (jWang et al.1 1 2 llh . The resulting confidence re- 
gions for all parameters arc not affected by thinning and 
it can be safely concluded that the steps in the MCMC 
runs are not subject to auto-correlation. All the accepted 
steps of each run are therefore used to create the poste- 
rior distributions. 

Following the first run, all assumptions (on distance, 
Nh and the presence of a PL) are progressively relaxed 



in the MCMC Runs #2 through #7, with the last one 
producing the i?NS measurement with the fewest assump- 
tions. The Mns - -Rns contours and i?NS distribution of 
each of these are displayed in Figures IMT51 and the re- 
sults are listed in Tables [5] and [6] In the process of re- 
laxing assumptions, one confirms that the results remain 
consistent between each run and such process does not 
significantly bias the results. The effects of relaxing each 
assumption are briefly described below. 

4.2.1. Using Gaussian Bayesian Priors for Distances 

When adding Gaussian Bayesian priors in place of the 
fixed distance parameters, the Mns-Rns contours are, 
as expected, broader in the i?oo direction. Because the 
normalization of a thermal spectrum such a nsatmos is oc 
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Fig. 12. — This figure is similar to Figure[9] but for Run #4. Here, the Gaussian Bayesian priors on the distances are included, as well 
as the PL spectral component, with Njj held fixed. The resulting NS radius is Rjqs = 8.0^j'?km. 



(Roc/d) , relaxing the assumptions on doc increases the 
possible values of Roo . This effect is mostly noticeable for 
the two targets observed with the highest S/N, i.e., M28 
and NGC 6397. In Run #2, the posterior distribution 
of i?NS corresponds to i?NS = 7.6l ; 9 km, broader than 
that of the previous run. 

4.2.2. Adding a Power-Law Spectral Component 

When adding PL components to account for possi- 
ble excesses of photons at high energy, one finds NS 
parameters consistent with those of the previous runs. 
Most PL normalizations are consistent with zero in 
Runs #3, #4 and #7. For M28, the PL normal- 
izations in these runs is consistent with zero, within 
2.4cr. For NGC 6397 in Runs #3 and #4, the con- 



sistency with zero is only marginal, within 2.8a. Fi- 
nally, in Run #7, the PL normalization of NGC 6397, 
NormpL.NGC 6397 = 2.7± 1 jxlO- 7 keV- 1 s- 1 cm- 2 , is 
not consistent with zero. This may indicate the possible 
contribution of a PL component at large energies and it 
is further discussed in Section [5j 



4.2.3. Relaxing the Nh assumption 

As the constraints on the Nh parameters are relaxed, 
the nsatmos best-fit parameters remain consistent with 
those of the previous runs. The posterior distributions 
of the five Nh parameters are consistent with the X-ray 
deduced Nh values found from the spectral fits of the 
sources individually, i.e., the values used in Run #1-4. 
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TABLE 6 

Results from the simultaneous spectral fitting, with free N h 



Target a pilcup kT cS M NS Roo N H ,22 PL Norm X 10 7 

(eV) (A/ ) (km) keV-^-^m" 2 

Run #5: Free N H , Fixed d,GC> No PL included, _R N s=7.5t ^'q km 

xg/dof (prob.) = 0.98/638 (0.66), 8% accept rate 

M28 0.43t£^ l-60t£^ 12.4+i;^ 0248^™^ 

NGC 6397 70t* 0.76±{J-^ 9.0±^ 0.105±g"^| 

M13 106^ l-43±°-;tL 11.4±?;? 0.010«;^ 

^Cen 137±1? 2.00±°-jf 16.7±1£ 0.152±°;°« 

NGC 6304 - 120±|| 1.13±g;|g p ljUlfj 0.315±g;ggg - 

Run #6: Free 2Vg , Gaussian Bayesian priors for dgc, No PL included, i?iMg=7.8^"j i km 

X^/dof (prob.) = 0.99/633 (0-59), 8% accept, rate 

M28 0l2™ 158^ i-SS-oisp 12 - 4 -!:? 0.248t£o£ 
NGC 6397 71±§ 0.78±°-:^ p 9-4±l"| 0.104±g gi| 

M13 l00l?| i.38±g;g p ll.3±|;| 0.0l0l°;^| p 

^Cen 133±|g 2.07±°;^ 17.11*;* 0.156±°;°*° 

NGC 6304 - U6±i4 l-09to.52 P 10 - 3 -i:9 °- 321 -o.06i - 

Run #7: Free , Gaussian Bayesian priors for doc ■ PL included, i?Ng=9.llj;'g km 
xl/dof (prob.) = 0.98/628 (0.64), 7% accept, rate 
0.34±g;it 137^1 UOT^L 12.6t^ O^is™ 5^ 
NGC 6397 67tg 0.86l°-* 7 10. 8l^ 7 0.1161°'°^ 2 -7±i "f 

M13 921?* 1.47tr4 12 - 6 -2'3 0.014lg;g» 2.4^ 

^Cen 130±g 2.42t°-«* 20.31^ 0.172l°^ 7 2.9j£», 

NGC 6304 112±?i l-32l°; 8 7 ? p 12.ll*;' 0.3461°;°^ ^-f.tp 



Note. — This tabic is similar to Tablc[5l abbreviations and symbols arc the same. The only difference 
is that the absorption Nh remained free for these runs. Quoted uncertainties arc 90% confidence as well, 
"p" indicates that the posterior distribution did not converge to zero probability within the hard limit of 
the model. 



TABLE 7 

Results from Simultaneous Spectral Fitting with nsagrav 

Target « pi i eup kT cS M NS Rco N H ,22 PL Norm xlO" 7 

(eV) (Mq) (km) keV" 1 s^ 1 cm" 2 

Run #8: Fixed N H , Fixed d G c , No PL included, i? NS =7.8lo 3 km 

xg/dof (prob ) = 1.04/643 (0.25), 10% accept, rate 

M28 0.46l£^ 166i^ 2 -03lo'.S P ;<i - ;! 'r," (0.252) 

NGC 6397 66±2 °- 51 -0A5 P s - 7 -o.t (0.096) 

M13 WltH 1.42±g;||^ 11.5t^ 7 2 (0.008) 

wCen 113^ 1.86+°;^ p 14.3jfj (0.182) 

NGC 6304 126±^ 1.4ll£?^ HA^ 1 (0.346) 



Note. — This table is similar to Table ]E\ But this run was performed with the nsagrav 
model instead of the nsatmos model for comparison purposes. The differences between the two 
models arc described in Scction l4.2l for example, in the high-rcdshift regime (Sec Figure fl6t . 
"p" indicates that the posterior distribution did not converge to zero probability within the 
hard limit of the model. Quoted uncertainties are 90% confidence. 

4.2.4. All Assumptions Relaxed 4.2.5. Comparison with nsagrav 



For this final run (#7), 7?ns = 9.izti"s ^ m ^ s consistent 
with the radii obtained in the previous MCMC runs (#1- 
6). The posterior distributions of A/ns _ ^?ns are shown 
in Figure [15] and detailed in Table |6] Once again, all 
resulting values are consistent with those of the previous 
runs. Progressively relaxing assumptions ensures a good 
understanding of the spectral fit, with no unexpected 
behavior. 



We also performed the fit using the model nsagrav 
instead of nsatmos for comparison purposes. It has pre- 
viously been shown that nsatmos and nsagrav produce 
simil ar spectral parameters when fit to experimental data 
(e.g. lWebb fc Barretl[2007l) . In XSPEC, for the nsagrav 
model, the R^s range is [6 km, 20 km] , and the Mns 
range is [0.3 Mq, 2.5 Mq] compared to [5 km, 30 km] and 
[0.5 Mq, 3.0 Mq] with nsatmos. This run (#8) was done 
with the same characteristics as Run #1. 
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Fig. 13. — This figure is similar to Figurc[9] showing results of Run #5. The hydrogen column density Njj is left free in this run, but 
the distances remained fixed, and no PL component was added. This run produced -Rns = km. 



When comparing the posterior distributions of the pa- 
rameters and A/ns~^?ns contours obtained with nsagrav 
(Figure [16] for Run #1 and Table [7]) to those obtained 
with nsatmos (Figure |9]and Table[5]), some consistencies 
can be noticed. However, not all distributions are consis- 
tent between the two models. Specifically, for M28, M13, 
and NGC 6304, one can notice that an additional distinct 
lobe at high Mns appears in the Mns~-Rns parameter 
space. This appears to be because the nsagrav model as 
implemented in XSPEC gives different values in this pa- 
rameter space than returned by nsatmos; the authors of 
this model state that this is because the model is inappli- 
cable in this parameter region 13 (Zavlin and Pavlov, priv. 

13 The A/ns-Rns space where 



comm.). For example, some sets of Mns _ ^?ns allowed 
by nsagrav and giving an acceptable x 2_vaiue lead to 
imaginary values of Roo. It is important for an observer 
to keep this fact in mind, otherwise, results produced 
by the XSPEC implementation of nsagrav could be mis- 
interpreted. In light of the pitfall mentioned here, the 
nsagrav model should be used with care. 

5. DISCUSSION 

This paper presented the simultaneous analysis of the 
spectra from five qLMXBs in GCs with a common i?Ns 
parameter for all targets. The posterior distributions 

nsagrav is applicable can be seen her e 

http: / /heasarc . gsf c .nasa. gov/xanadu/xspec/models/m-r .pdf 



24 




Rn$ {km) 



K NS (km) 



Fig. 14. — This figure is similar to Figure [9] but for the MCMC Run #6. The characteristics of this run include Nh values free to 
vary in the fit, and the presence of the Gaussian Bayesian priors for the distances. No PL component was included in run. A value 
-Rns = 7.8^-,'? km was found. 



for i?NS> -^NSj Roo, kT e g, and Nh were obtained from 
MCMC simulations, which included Gaussian Bayesian 
priors for the distances to the GCs hosting the targets. In 
this discussion section, the original work performed here 
and the data used are summarized. This is followed by 
a subsection discussing various possible biases resulting 
from the MCMC analysis. The discussion finishes with 
the implication that the resulting i?NS measurement may 
have for the determination of the dEoS. 

5.1. List of New Analysis Methods, Data and Results 

The following two paragraphs aim at summarizing the 
novel approach to the analysis of the NS thermal spec- 
tra, and unused data presented in this paper. The 



MCMC framework for spectral analysis is a rather re- 
cent approach, made convenient by the development of 
PyXSPEC, the Python version of XSPEC. The "Stretch- 
Move" MCMC algorithm used here differs from the usual 
Metropolis- Hasting (M-H) algorithm of standard MCMC 
simulations. It has been developed recently and while 
it presents significant advantages over M-H, it is still 
scarcely used in astrophysics. This work includes X-ray 
data not previously presented in the literature, namely 
the Chandra exposure of NGC 6304. 

This paper also contains a more complete presentation 
of the resulting best- fit NSs physical parameters (kT c g, 
Mns, and i?Ns)- Specifically, the products of the MCMC 
simulations in 2-dimensional matrices of posterior dis- 
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Fig. 15. — This figure is similar to Figure[9] corresponding to Run #7. Here, all the possible assumptions have been relaxed to obtain a 
-Rns measurement the least affected by systematic uncertainties. The Nh parameters are left free; and Gaussian Bayesian priors and PL 
components are included. This results in an i?NS measurement: i?NS = 9.1^'gkm 



tributions (see Figures HH5]) are displayed. Such figures 
better represent the true distributions of the NS physi- 
cal parameters compared to simple lists of best-fit values 
with their uncertainties, mostly because the distribution 
are not necessarily Gaussian (particularly with models 
like nsatmos). We encourage researchers to present their 
results in such a way. 

Finally, t he approach of this paper i s not dissimilar to 
the work of ISteiner et al.l (|2010t 120121 ) in the sense that 
different targets are combined to produce constraints 
on the dEoS. However, our analysis imposes i?NS to be 
quasi-constant for all the targets, as justified by recent 
observations favoring "normal matter" dEoSs. Other 
differences include two qLMXBs (M28 and NGC 6304) 



added to the present analysis, and the qLMXB in 47 Tuc 
not used here because of the uncertaint i es related to pile- 
up. Finally, the work of ISteiner etaD (|2010l 12012( 1 uses 
type-I X-ray burst sources, which are not considered in 
the present analysis. 

5.2. Possible biases resulting from the analysis 
5.2.1. Non-zero power-law for NGC 6397 

The runs that included a PL component (#3, #4, 
and ff7, with photon index r = 1) resulted in PL 
normalizations consistent with zero for all targets, ex- 
cept for NGC 6397. Specifically, NormpL,NGC 6397 = 
2.7ti'3><10 -7 keV -1 s _1 cm~ 2 , which is not consistent 
with zero, at the 3.4c level, in Run ffl . This in- 
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Fig. 16. — This figure is similar to Figure[9] i.e., Run #1 with fixed distances instead of Gaussian Bayesian priors, fixed Nh values, and 
no PL component included. The only difference resides in the H-atmosphere model, where nsagrav model was used instead of nsatmos. 

is consistent with that of Run #1, but 
This is further discussed in Section l4.2l 



Note that the lower limit for -Rns m nsagrav is 6 km. The i?NS measurement, -Rns = 7.8_ ' 3 km, 
the some of the Mns~-Rns posterior distributions are significantly different from those of Figure[2] 



dicates a possible non-negligible contribution of a PL 
above - 2keV. Using XSPEC (without the MCMC 
approach), the spectra of NGC 6397 are fitted without 
the other qLMXBs, and the PL normalization obtained 



Nonetheless, adding a PL contribution for each tar- 
get does not significantly bias the posterior distribu- 
tion of i?Ns- Specifically, adding PL components be- 
tween Runs #1 and #3 changed the radius measurement 



is Normp L = 2.1^3 xlO-^keV" 1 s" 1 cm -2 which corre- from i? NS = 7.1±°;g km (Run #1) to R 



-1.3 

sponds to a contribution of 2.61^% of the total unab- 
sorbed flux of NGC 6397 in the 0.5-10 keV energy band. 
Such a contribution is consiste nt with that measur ed in 
the previous results, < 3.3% (jGuillot et all 1201 lal) . for 
the same photon index. 



NS 



7.3t° ; 5 6 km 



(Run #3); between Runs #2 and #4, JSns changed from 
i?NS = 7.6±o:9km (Run #2) to i? N s = 8.0l^km (Run 
#4). Between runs #6 and #7 (free Nh), the i?NS 
distributions changed more significantly, but they are 
nonetheless consistent: -Rns = 7.81^ km (Run #6) to 
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R NS = 9.lt\ l km (Run #7). Therefore, adding PL com- 
ponents does not significantly change the i?NS posterior 
distribution, but nonetheless includes systematic uncer- 
tainties related to the possible presence of a PL compo- 
nent into the measured i?NS- With the limited PL con- 
tributions observed (< 5%), the choice of photon index 
does not affect the nsatmos component. This is tested 
on NGC 6397 alone, the qLMXB with the strongest PL 
contribution, where the PL photon index is changed from 
L = 1 to L = 2. Such a change results in consistent PL 
contributions, but more importantly, it results in a non- 
significant increase of only ~ 1% in the value of i?NSi 
well within the uncertainties of the measurement. Other 
parameters such as the mass, temperature and the PL 
normalization are consistent between the two trials, with 
L = f and with L = 2. Given that NGC 6397 has the 
strongest PL contribution and that changing the photon 
index from L = 1 to L = 2 does not modify -RnSj it is 
expected to be the same for the simultaneous spectral 
analysis. 

5.2.2. Effects of Individual Targets on the Simultaneous Fit 

An additional set of MCMC runs was performed to 
investigate whether some targets have a dominant bias- 
ing effect on the simultaneous spectral fitting presented 
above. This analysis is performed by excluding each tar- 
get individually and observing the resulting marginalized 
posterior distributions of the parameters. Results are 
listed in Table [HI These runs were performed like Run 
#7, i.e., with the Nh value free, with Gaussian Bayesian 
priors on the distances, and with the additional PL com- 
ponent. The values of i?NS obtained when removing each 
target are all consistent within 2er of each other, and 
more importantly, consistent with the i?NS distribution 
obtained when all targets are included (Run #7, in Ta- 
ble [S]). Similarly, the kT e s, Mns, Roo and the a p ii eU p val- 
ues of each targets are consistent with each other in all 
five cases presented in Table HI and also consistent with 
the values from Run #7. Overall, no individual target 
has any dominant effect on the simultaneous fit. This en- 
sures that no qLMXB induces a significant bias on i?NS- 
This is also compatible with recent results (jSteiner et al.l 
12012ft . obtained by combining Mns~-Rns distributions of 
qLMXBs and type-I X-ray bursts. In that work, it was 
demonstrated that removing extreme cases (for example 
the qLMXB X7 in 47Tuc, or the qLMXB in M13) had 
no or little effect on the resulting empirical dEoS. 

For completeness, and to confirm this observation, a 
spectral analysis combining X7 to the other five qLMXBs 
of this work is performed in XSPEC in order to deter- 
mine i?NS- Note that the MCMC approach is not used 
and is not necessary here since this additional analysis 
is simply a consistency check to determine if outliers 
can have an impact of the best-fit i?NS m easurement 14 . 
Furthermore, Gaussian Bayesian priors are not used for 
this consistency check. The spectral data available for 
X7 and used here are described in a previous work 
(|Heinke et all 120061) . For the distance to the qLMXB, 

14 As mentioned before, X7 was not used in the main analysis 
of the present work because of the unquantified systematic errors 
associated with the correction of the large amount of pile-up af- 
fecting the spectra, and its uncertain effect on our final statistical 
error bars. 



TABLE 8 

Effect of Individual Targets on the Simultaneous Spectral 

Fit 



Target 






xS/d-o.f. (prob.) 


Accept. 


excluded 


km 












rate 


NONE (Run #7) 


9 1+ 1 - 3 

a - -1.5 


km 


0.98 


/ 


628 


(0.64) 


7% 


WITHOUT M28 




km 


0.98 


/ 


381 


(0.69) 


7% 


WITHOUT NGC 6397 


10.7th 7 


km 


0.89 


/ 


428 


(0.95) 


9% 


WITHOUT M13 




km 


0.94 


/ 


588 


(0.86) 


7% 


WITHOUT u Cen 




km 


0.95 


/ 


601 


(0.81) 


8% 


WITHOUT NGC 6304 


9.0±l : g 


km 


0.93 


/ 


622 


(0.88) 


8% 



Note. — The spectral fits in this table were performed following 
Run #7, with free Nh values, Gaussian Bayesian priors on the distances 
and with an additional PL components in the model. Each target were 
successively removed to investigate the possible effect of individual sources 
on the global fit. The i?NS values obtained in each of these five tests are 
consistent with each other, and with that of Run #7. This confirms 
that none of the five qLMXBs significantly skews the -Rns measurement. 
Quoted uncertainties arc 90% confidence, and values in parentheses are 
held fixed. 

we used the weighted average of a ll recent distance est i- 
mates listed in another reference (jWoodlev et al.i r2012). 
i.e., d47Tuc = 4.52 kpc. The resulting best-fit NS radius 
found was i?NS = ^-^-2 0p ^ m with the acceptable statis- 
tics xl/doi (prob.) = 0.98/1003 (0.70). This result is 
consistent with the i?NS measurement of Run #1 (per- 
formed without adding the spectra of X7) and seem to 
confirm that outliers, such as the qLMXB X7 in 47Tuc, 
do not affect the r adius measurement, as demonstrated 
in a previous work (jSteiner et alj|2012l ). 

5.2.3. Composition of the Neutron Star Atmosphere 

The Mns - -Rns measurements of NSs in qLMXBs rely 
on the atmosphere modeling, which itself relies on a ma- 
jor assumption of this work, namely the composition 
of the NS atmosphere. It is generally assumed that 
NS atmo sphere are composed exclusively of pure hydro- 
gen (e.g.. lRutledge et al.l[200l iHeinke et al.ll2006f ). since 
heavier accreted elements will settle through th e atmo- 
sphere on short time scales (Bildsten et al.lll992l ). NS He 
atmospheres could be observed in the case of He accretion 
from WD donors in ultra-compact X-ray binaries. An 
He-atmosphere model was used to fit the spectra of M28 
and led to a radius 1.5 x larger than that obtained with 
a H- atmosphere, favoring stiff dEoSs (Servillat et al. 
12012ft . Similarly, for the qLMXB in M13, assuming a 
He composition of the NS atmosphere increases th e ra- 
dius i?NS by a factor ~ 1.2 l)Catuneanu et al.ll2013l) . 

The actual composition of the NS atmosphere can be 
infered from the identification of the donor companion 
star, which proves a difficult task in the crowded envi- 
ronment of GCs. Only two GC qLMXBs have iden tified 
counterpa rts, X5 in 47Tuc ([Edmonds et all 12002). and 
in to Cen ( Haggar d et al.l I2004D , both of which discov- 
ered from their strong Ha emission, indicating hydrogen 
donor stars. 

5.2.4. Causality Limit 

With the assumptions made and the model chosen in 
the present work, parts of the Mns-Rns contours result- 
ing from the analysis cover a section of the parameter 
space that goes past the causality limit as set in ealier 
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works (jLattimer et al.lll990t ILattimer fc Prakashll200lD . 

Stricter constraints on the Mns~Rns contours could be 
obtained by imposing that the sets of MCMC accepted 
points (-Rns^ns) all obey i?NS > 3GMns/c 2 , i.e., do 
not cross the causality line. However, we choose not to 
change constraints; this produces a larger error region 
than may be necessary (if one were to explicitly adopt 
causality as an assumption) , but the goal of this analysis 
is to produce the most conservative, assumption-free un- 
certainty region for i?NS- Only if much higher S/N data 
were obtained, and the M-R parameter space required (or 
strongly preferred) a value in the r egion excluded by the 
causality requirement cited above (jLattimer et al.lfl990f) 
would it become necessary to revisit this assumption. 

5.2.5. Effect of Assumptions 

As discussed in this paper, assumptions can have a 
strong effect on the interpretation of spectral fits of 
individual sources and consequently on the simultane- 
ous spectral analysis as well. The selection of the dis- 
tances to the GCs (fixed or with Gaussian Bayesian pri- 
ors) can skew the i?NS measurement toward smaller or 
larger values. For instance, in early runs of this analy- 
sis, th e distance du2s — 5.1±0.5kpc (H ees fc Cudworthl 
119911 ) was initially used, before it was updated to a more 
recent value dM28 = 5.5 ± 0.3 kpc (|Testa et all [20011: 
iServillat et all 12012ft . This caused the i? N s measure- 
ment (of Run #7) to change from i?NS = 8.7^^ km 
to i?NS = 9.1^}'5km as the distance of M28 was in- 
creased. The values are consistent with each other, but 
this larger d,M2s shifted the i?NS measurement. This is 
expected since the normalization of a thermal spectrum 
is oc (Roo/df. 

Another strongly influential parameter is the galactic 
absorption. A change in Nh will strongly affect the Mns 
and J?ns best-fit. Specifically, a decrease (increase) in 
Nh results in an decrease (increase) in R x , respecti vely. 
This was shown for NGC 6397 and cj Cen, in Section |4~T1 
For NGC 6397, the Hi-survey N H value (N H ,22 = 0.14) 
leads to a Roa — 11. 8^7 km, while the X-ray deduced 
value {N H ,22 = 0.096±g;g^) produces i?^ = 8.4±^km, 
a decrease of ~ 30%. For oj Cen, the best- fit Roo al- 
most doubles as the constraint on Nh is relaxed, causing 
the value to increase from Nh,22 = 0.09 (Hi-survey) to 
Nh,22 = 0.182+q;q47. Consequently, using assumptions 
for the Nh values, like those derived from Hi-surveys, 
can lead to strongly skewed i?NS measurements. 

Consequently, it is preferable to avoid using assumed 
values, when possible, and therefore the i?NS measure- 
ment from Run #7 is presented as the final result of this 
work. 

5.3. The i?NS measurement 

The goal of this paper being to measure -Rns, the fol- 
lowing discussion is focused on the i?NS posterior dis- 
tributions, on the comparison with other i?NS measure- 
ment, and on the implication that such a radius mea- 
surement will have for the determination of the dEoS. 

The striking observation one can make pertains to the 
low range of values of the resulting _Rns distributions 
obtained from the different runs. The .Rns distribu- 
tions remained below R^s < 10.4 km (90%-confidence) , 



or i?NS < H.lkm (99%-confidence). This resulted from 
Run #7 (see Figure IT5|). where a particular effort was 
made to consider all possible sources of systematic uncer- 
tainties. The result from Run #7 is the most general .Rns 
distribution, i.e., with the fewest assumptions, that can 
be produced. Also, the progressive relaxation of the as- 
sumptions throughout the analysis demonstrated that no 
unexpected behavior was present in the final Mns~Rns 
distributions of Run #7 and that the resulting low-value 
of i?NS was not affected by systematics. 

Previous works reported low values of NS radii, but 
these measurements have high uncertainties due to low 
S/N, leadin g to poorly constrained i?NS and Mns ( e -R-i i n 
NGC 2808. IWebb fc Barred [20071 IServillat et all 12008ft . 
Another qLMXB in NGC 6553 was identified with a 
small radius, -Rn s = 6.3tn'sk m (90%-confidence) for 
Mns = 1.4 M (|Guillot et all 1201 lbft . However, low- 
S /N Chandra observations demonstrated that the XMM 
spectra of the source was affected by hard X-ray contami- 
nation from a marginally resolved nearby source. Higher- 
S /N observations with Chandra are necessary to confirm 
the qLMXB classification and produce the uncontami- 
nated spectrum necessary for its use in the present anal- 
ysis. 

In addition to qLMXB _Rns measurements, low radii 
were found from the analysis of photospheric radius ex- 
pansion type-I X-ray bursts. A review of the method 
used to determine i?NS from these sources can be found 
in the literature (jOzell 120061: iSuleimanov et all l2011bt) . 
The LMXBs EXO 1745-248, 4U 1608-52, and 4U 1820- 
30 were found to have respective radii in the 2a ranges 
i?NS = [7 5 - 11.0] km (lOzel et al.| [2009l). i? N s = [7.5 - 
11.51 km dGuyer et al.l l2010aft and R NS = [8.5 - 9.5] km 
(jGiiver et al.l l2010bf ). respectively. While these results 
are on a par with what is found in this paper, controversy 
emerged with the realization that the analysis presented 
in the cited works was not internally consistent because 
the most probable observables (from Mon te- Carlo sam- 
pling ) led to imaginary masses and radii ([Steiner et al.l 
2010;. Relaxing the assumption that the photospheric 
radius equals the physical radius i?NS at touchdown led 
to real- valued solutions of Mns and i?NS> and to larger 
upper limits for the radius. Furthermore, it is argued 
in a later work that the short bursts from EXO 1745- 
248, 4U 1608-52 and 4U 1820-30 are not appropriate for 
such analysis because the post-burst cooling evolution 
of these sour ces does not match the th eory of passively 
cooling NSs (jSuleimanov et all l2011af ). Therefore, the 
Mns _ -Rns constraints from type I X-ray bursts should 
be considered with these results in mind. 

More recently, distance independent constraints in 
Mns~-Rns space were produced from the analysis 
of the sub-Eddington X-ra y bursts from the type I 
X-ray burster GS 1826-24 (jZamfir et al.ll20Tl . That 
analysis, performed for a range of surface gravities 
(log 10 (g) = 14.0, 14.3, 14.6) and a range of H/He abun- 
dances (0.01 Z Q , 0.1 Z Q and Z Q ) led to radii i?NS < 
11.5 km. While distance-independent, the results are 
highly influenced by the atmosphere composition and 
metallicity. For pure He comp osition, the upper limit 
of i?NS becomes i?NS < 15.5 km (jZamfir et al.N2012l ). 

Finally, the multiwavelength spectral energy distri- 
bution of the isolated neutron star RX J 185635-3754 



The Radius of Neutron Stars 



29 




tf N s (km) 



+ 1.3 



km 



Fig. 17. — This figure shows the constraint on the dEoS imposed by the radius measurement obtained in this work: i?NS = 9.1_ 1 5 
(90%-confidence). The dark and light shaded areas show the 90%-confidence and 99%-confi dence constraints o f the -Rns measurement, 
respectively. The mass measurement of PSR J1614-2230 is shown as the horizontal band (Demorcst ct al. 2010). "Normal matter" EoSs 
are the colored solid lines. Other types of EoSs, such as the hybrid or quark-matter EoSs are included for comparison, with dashed lines. 
As mentioned in Section [3] the present analysis only places constraints on the "normal matter" EoSs since t hey are the only family of EoSs 
included in our assumptions. Among the m, only the very soft dEoSs (such as WFF1, IWiringa et al.lH98SI) are consistent with the radius 
obtained here. The EoS are obtained from Lattimcr & Prakash (2001, 2007). 



was analyzed to produce small values of i?NS and -Mns 
with no plausible dEoS consistent with these values: 
i?NS ~ 6fe m a . n d -Wns ~ 0.9 M for d = 61 pc 
(|Pons et al.l 120021) . A recent distance estimation to the 
source d = 1231^ pc (jWalter et al.ll2010D led to revised 
values: R NS = 11 .5±1.2km and M NS = 1.7±1.3M Q 
(jSteiner et al.l l2012'). While this result is consistent with 
the i?NS measurement obtained in this paper and with 
the other works reporting low-i?NS values, it has to 
be taken with care since the high-magnetic field of the 
source is not accounted for in the spectral model used by 
the original analysis. 

Recently, it was shown that the dEoS can be empir- 
ically determined from Mns--Rns measurements of NS, 
using the thermal spectra of qLMXBs a nd the photo- 
spher ic radius expansion of X-ray bursts (jSteiner et al.l 
I2010D . This method uses MCMC simulation and 
Bayesian priors to determine the most probable dEoS 
parameters, and equivalently, the corresponding most 
probable Mns{Rns) for NS. In a recent paper, this 
method was used with four X-ray bursting sources and 
four GC qLMXBs. Considering all scenarios, the 2a 
lower and upper lim its for i?NS are 9.17 km and 13.92 km 
(|Steiner et al.l l2012). The i?NS distribution of the present 
paper i?NS = 9.1+^5 km (90%-confidence, from Run 
#7) is consistent w ith several of the model variations of 
iSteiner eT al. (2012), namely variation C (dEoS param- 
eterized with uniform prior in the pressure at four en- 
ergy density values), variation CII (same as previous, but 
with low value of the color correction, 1 < fc < 1.35), 
variation AII/AIII (dEoS parameterized as two piece- 
wise continuous pow er-laws, with 1 < /c < 1-35), see 
ISteiner et al.1 (|2012f ) for details about the variations of 
the model. Variation E (dEoS for quark stars) is incom- 



patible with our original assumption that i?NS is quasi- 
constant for a large range of Mns above 0.5 M©. 

Theoretical EoSs have been proposed for more than 
two decad es. A non-exhaustive list can be found in the 
literature (jLattimer fc Prakash|[200lll2007D . When com- 
paring the resulting i?NS distribution to proposed theo- 
retical "normal matter" dEoSs, one can note that most 
of those are not consistent with the low-i?NS result pre- 
sented in this work. Indeed, most of the dEoSs describing 
"normal matter" correspond to radii larger than 11.5 km 
(see Figure[17]). A spread in i?NS ia observed in these dE- 
oSs at large masses, in the part of the Mns _ ^?ns diagram 
where the compact object approaches collapse. However, 
this breadth of the i?NS variation for a given dEoS is well 
within the uncertainties obtained in this work. Overall, 
the radius measurement i?NS = 9.lj^'g km constrains the 
dEoS to those consi stent with low-i?NS, such as WFF1 
(jWiringa et al.l[T988| ) . Note that this analysis cannot ad- 
dress the veracity of more exotic types of EoSs (hybrid 
and SQM) or any dEoS which does not predict a quasi- 
constant i?NS within the observable mass range. 

It is known that i?NS is related to fundamental nu- 
clear physics parameters, such as the symmetry en- 
ergy (Horowit z fc Piekarewiczll2001bl lal). We expect the 
present constraints on i?NS can be used to constrain this, 
and other properties of dense nuclear matter. We leave 
this for future work. 

It has been pointed out ([Lattimer fc Prakash! I2010D 
that an a rgument regarding a maximally compact neu- 
tron star (jKoranda et al.lfl997| ) results in a relationship 
between the maximal neutron star radius (i? max ) and the 
maximal neutron star mass (M max ) for a given equation 
of state: 



Mrr 



G 

2.824— 



(4) 
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where G is Newton's constant and c is the speed of light. 
Adopting this, the 99% confidence upper limit is M max < 
2.66 Mq, which does not violate any measured neutron 
star masses at present. 

The small- i?NS value found in this paper, and other 
low-i?NS results cited above, are consistent with soft 
dEoSs such as WFF1. However, results of Mns and 
i?NS measurements from other sources seem to fa- 
vor stiffer dEoSs. The qLMXB 47Tuc X7 has a re- 
ported radius i? N s = 14.5j^;|km for Mns = 1.4 M© 
(iHeinke et al.l 12006ft. suppor ting stiff dEoSs, such as 
MSO/2 (jMuller fc Serotlll996ft . Nonetheless, the range of 
radii allowed by the published Mns _ -Rns contours for X7 
is consistent with the radius measurement presented in 
the present work. Moreover, the X7 Mns ~-Bnr contours 
are c ompatible with the dEoS WFF1 ( Wiringa et al.l 
1988). Another work used the long photospheric radius 
expansion X-ray bursts from 4U 1724-307 to conclude 
that stiff dEoSs are describin g the dense matter inside 
NSs (|Suleimanov et al.ll2011al ). Specifically, it was found 
that i?NS > 13.5 km for Mns < 2.3 M , and for a range 
for NS atmospheric composition. Lower _Rns values, in 
the range 10.5-17 km, are allowed for Mns > 2.3 M Q , 
for pure H or solar metallicity composition. This ra- 
dius measurement is only marginally consistent with the 
present work for large masses, Mns > 2.3 M Q , which 
implies a dEoS capable of reaching A/ns ~ 2.3 M for 
i?NS ~ 10 — 11km. Finally, another radius measure- 
ment, obtained by modeling the ther mal pulses of th e 
millisecond pulsar PSR J0437-4715 (|Bogdanovi [20121 ) . 
led to values, i?NS > H km (3a), is inconsistent with the 
measurement presented in our work. 

6. SUMMARY 

In this paper, we measured i?NS using the assump- 
tion that the radius is quasi-constant for a wide range 
of Mns larger than Mns > 0.5 M Q , i.e., constant within 
the measurement precision. This is justified by recent 
observations favoring "normal matter" dEoSs which are 
described by this characteristic. For this analysis, the 
spectra from five GCs qLMXBs observed with the Chan- 
dra X-ray Observatory and XMM-Newton were used in a 
simultaneous analysis, constraining i?NS to be the same 
for all targets. 

For this, we used an MCMC approach to spectral fit- 
ting, which offers several advantages over the Levenberg- 
Marquardt ^-minimization technique generally used for 
spectral fits. For example, the MCMC framework al- 
lows imposing Bayesian priors to parameters, namely the 
distance to the host GCs. By doing do, the distance 
uncertainties are included into the posterior i?NS distri- 
bution. In addition, one can marginalize the posterior 
distributions over any parameters and very easily obtain 
Mns _ -Rns distributions, while the grid-search method in 
XSPEC can be problematic in the case of spectral fits 
with many free parameters and complicated x 2 -space. 



The algorithm chosen in this work is an affine-invariant 
ensemble sampler, commonly called "Stretch-Move" al- 
gorithm, which is particularly appropriate (i.e., converg- 
ing efficiently) for elongated and curved distributions. 

The principal result of the simulations performed in 
this analysis is that NSs are characterized by small phys- 
ical radii. Specifically, when the distances and Galac- 
tic absorption parameters are fixed, i?NS = 7.1^Q gkm 
(from Run #1). A more general posterior distribution 
for _Rns i.e.) less prone to systematic biases, is obtained 
by applying Gaussian Bayesian priors for the five GC 
distance, by freeing the Nh parameters, and by adding 
a PL component to the model to account for a possi- 
ble spectral component at high photon energies. Such a 
spectral component could be the largest possible source 
of uncertainty, and could be skewing _Rns downward, but 
it is accounted for in the last and most relaxed MCMC 
run. In fact, such a spectral component was discovered 
herein for NGC 6397. 

The progressive relaxation of assumptions led us to a 
good understanding of the spectral fit in Run #7, min- 
imizing systematic uncertainties. Therefore, with the 
H-atmosphere model nsatmos, the measured NS radius 
is i? N s = 9.ltijkm (from Run #7). These results 
are compatible with other low-i?NS measurements from 
GC qLMXBs or type-I X-ray bursts, but not consistent 
with some published i?NS measurement leading to values 
-Rns > Hkm. We recommend these -Rns constraints, 
from Run #7, be those relied upon for constraints on 
the dEoS and other nuclear physics model parameters, 
as this run has the fewest associated assumptions behind 
it. 

Among the dEo S listed in previous works 

(jLattimer fc Prakashl I2001L l2007t) . the i? N s mea- 
surement presented here is only compatible with 
"normal matt er" dEoSs consistent with i?NS ~ 10 km, 
e.g. WFF1 (jWiringa et al.l 119881 ). Most dEoSs are 
compatible with larger radii, at i?NS ~ 12 km and above. 
Given the results presented in this work, the theory of 
dense nuclear matter may need to be revisited. 
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